Optimization with LBFGS gives DimensionMismatch("dimensions must match")

Specific Domains Optimization (Mathematical)
,
1

My objective function has to solve an ODE. When I use the gradient-free algorithm (e.g NelderMead) to do optimization, everything is OK. But with gradient-based (e.g. LBFGS) I get the error: DimensionMismatch(“dimensions must match”)

It happens only if a parameter is a vector, in case of scalar I don’t have an error.

I have a small example to illustrate an issue.

using DifferentialEquations, Optim, LinearAlgebra

function f(dxdt,x,p,t)
    dxdt[1] = -p[1]*x[1] + p[2]
    dxdt[2] = p[1]*x[1]
end

tspan = (0.,8.)
tdata = range(0.,stop=8.,length=10)
p0 = [1., 1.0]
x0 = [10.,0.]

prob = ODEProblem(f,x0,tspan,p0)

sol = DifferentialEquations.solve(prob; saveat=tdata)

xdata = sol[2,:] + randn(length(sol[2,:])) # in-silico data

function loss_function(x)
    _prob = remake(prob; u0=convert.(eltype(x), prob.u0), p=x)
    sol_new = DifferentialEquations.solve(_prob; saveat=tdata)
    norm(sol_new[2,:] - xdata)^2
end

optimize(loss_function, [2.0, 3.0], NelderMead())    # OK
optimize(loss_function, [2.0, 3.0], LBFGS()) # ERROR

Error:

DimensionMismatch("dimensions must match")

Stacktrace:
 [1] promote_shape at .\indices.jl:154 [inlined]
 [2] promote_shape at .\indices.jl:145 [inlined]
 [3] -(::Array{Float64,1}, ::Array{Float64,1}) at .\arraymath.jl:38
 [4] loss_function(::Array{Float64,1}) at .\In[21]:15
 [5] finite_difference_gradient!(::Array{Float64,1}, ::typeof(loss_function), ::Array{Float64,1}, ::DiffEqDiffTools.GradientCache{Nothing,Nothing,Nothing,Val{:central},Float64,Val{true}}) at C:\Users\mzhen\.julia\packages\DiffEqDiffTools\bzjuh\src\gradients.jl:282
 [6] (::getfield(NLSolversBase, Symbol("#g!#15")){typeof(loss_function),DiffEqDiffTools.GradientCache{Nothing,Nothing,Nothing,Val{:central},Float64,Val{true}}})(::Array{Float64,1}, ::Array{Float64,1}) at C:\Users\mzhen\.julia\packages\NLSolversBase\KG9Ie\src\objective_types\oncedifferentiable.jl:56
 [7] (::getfield(NLSolversBase, Symbol("#fg!#16")){typeof(loss_function)})(::Array{Float64,1}, ::Array{Float64,1}) at C:\Users\mzhen\.julia\packages\NLSolversBase\KG9Ie\src\objective_types\oncedifferentiable.jl:60
 [8] value_gradient!!(::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}) at C:\Users\mzhen\.julia\packages\NLSolversBase\KG9Ie\src\interface.jl:82
 [9] value_gradient!(::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}) at C:\Users\mzhen\.julia\packages\NLSolversBase\KG9Ie\src\interface.jl:69
 [10] value_gradient!(::Optim.ManifoldObjective{OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}}, ::Array{Float64,1}) at C:\Users\mzhen\.julia\packages\Optim\0je8H\src\Manifolds.jl:50
 [11] (::getfield(LineSearches, Symbol("#ϕdϕ#6")){Optim.ManifoldObjective{OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}},Array{Float64,1},Array{Float64,1},Array{Float64,1}})(::Float64) at C:\Users\mzhen\.julia\packages\LineSearches\WrsMD\src\LineSearches.jl:84
 [12] (::LineSearches.HagerZhang{Float64,Base.RefValue{Bool}})(::Function, ::getfield(LineSearches, Symbol("#ϕdϕ#6")){Optim.ManifoldObjective{OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}},Array{Float64,1},Array{Float64,1},Array{Float64,1}}, ::Float64, ::Float64, ::Float64) at C:\Users\mzhen\.julia\packages\LineSearches\WrsMD\src\hagerzhang.jl:140
 [13] HagerZhang at C:\Users\mzhen\.julia\packages\LineSearches\WrsMD\src\hagerzhang.jl:101 [inlined]
 [14] perform_linesearch!(::Optim.LBFGSState{Array{Float64,1},Array{Array{Float64,1},1},Array{Array{Float64,1},1},Float64,Array{Float64,1}}, ::LBFGS{Nothing,LineSearches.InitialStatic{Float64},LineSearches.HagerZhang{Float64,Base.RefValue{Bool}},getfield(Optim, Symbol("##22#24"))}, ::Optim.ManifoldObjective{OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}}) at C:\Users\mzhen\.julia\packages\Optim\0je8H\src\utilities\perform_linesearch.jl:40
 [15] update_state!(::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Optim.LBFGSState{Array{Float64,1},Array{Array{Float64,1},1},Array{Array{Float64,1},1},Float64,Array{Float64,1}}, ::LBFGS{Nothing,LineSearches.InitialStatic{Float64},LineSearches.HagerZhang{Float64,Base.RefValue{Bool}},getfield(Optim, Symbol("##22#24"))}) at C:\Users\mzhen\.julia\packages\Optim\0je8H\src\multivariate\solvers\first_order\l_bfgs.jl:198
 [16] optimize(::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}, ::LBFGS{Nothing,LineSearches.InitialStatic{Float64},LineSearches.HagerZhang{Float64,Base.RefValue{Bool}},getfield(Optim, Symbol("##22#24"))}, ::Optim.Options{Float64,Nothing}, ::Optim.LBFGSState{Array{Float64,1},Array{Array{Float64,1},1},Array{Array{Float64,1},1},Float64,Array{Float64,1}}) at C:\Users\mzhen\.julia\packages\Optim\0je8H\src\multivariate\optimize\optimize.jl:57
 [17] optimize(::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}, ::LBFGS{Nothing,LineSearches.InitialStatic{Float64},LineSearches.HagerZhang{Float64,Base.RefValue{Bool}},getfield(Optim, Symbol("##22#24"))}, ::Optim.Options{Float64,Nothing}) at C:\Users\mzhen\.julia\packages\Optim\0je8H\src\multivariate\optimize\optimize.jl:33
 [18] #optimize#87 at C:\Users\mzhen\.julia\packages\Optim\0je8H\src\multivariate\optimize\interface.jl:116 [inlined]
 [19] optimize(::Function, ::Array{Float64,1}, ::LBFGS{Nothing,LineSearches.InitialStatic{Float64},LineSearches.HagerZhang{Float64,Base.RefValue{Bool}},getfield(Optim, Symbol("##22#24"))}, ::Optim.Options{Float64,Nothing}) at C:\Users\mzhen\.julia\packages\Optim\0je8H\src\multivariate\optimize\interface.jl:115 (repeats 2 times)
 [20] top-level scope at In[24]:1
2

I’m precompiling DifferentialEquations :slight_smile: but did you actually check the size of sol[2,:] and xdata in your loss call? That’s what it’s saying has incompatible dimensions.

3

Thank, you @pkofod .
You were right.

It started from correct size, but after few iterations size(sol_new[2,:]) becoming 7 instead of 10!!!

Does it mean that ODE solver cannot solve equations and just skips some time-points?

4

Sorry, I don’t know much about the ode solution interface, we’ll have to wait for someone else :slight_smile:

5

I encounter the same problem with BFGS, have you found a solution back then @mzhenirovskyy ?

6

This has nothing to do with BFGS. The issue is that the loss function was not written to be robust to divergent trajectories. It’s missing a check for sol_new.retcode == :Success. See Handling Divergent and Unstable Trajectories · SciMLSensitivity.jl

1 Like
7

I have the same problem, has this somehow been addressed?

8

Given that this thread is solved not by anything in the library but was a user code issue, yes it was addressed before the thread was opened in 2019 using the return codes. Do you have the same issu and did you try the same solution?

9

Hi @ChrisRackauckas Thank you for answering!

Yes, I do have the same issue. The website you linked to throws a 403 error. Is there a description of how to fix that? When I use Newton’s method it works fine by the way.

10