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

differentiation
ode
#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
0 Likes

#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.

0 Likes

#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?

0 Likes

#4

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

0 Likes