Using StaticArrays in NewtonMethod() optimization with simultaneous cost, gradient and hessian calculation

[fb456]

Hi,

I'm using Julia Version 1.8.5 and Optim v1.7.4
I am trying to optimize the following function, where a second order method is employed and the cost, gradient and hessian are calculated simultaneously. x is an array of length 2.
optimize(only_fgh!((F,G,H,x) -> costgradhess(F,G,H,x,x_0,y_0,a,b,alpha)), x_in, Newton())

where x_0,y_0,a,b,alpha are parameters needed by costgradhess. The optimization process worked perfectly when I use Vector{Float64} type for x. However, I decided to use StaticArrays (MArrays to be specific) for x to increase performance. Here is what it looks like (the code is quite long, so I skipped some parts)

function costgradhess(F,G,H,x_in::MVector{2,Float64},x_0::SVector{3, Float64},y_0::SVector{3, Float64},
a::SVector{3, Float64},b::SVector{3, Float64},alpha::SVector{3, Float64})

    if !(isnothing(G))
         # G[1] = ... ; G[2] = ...
    end
    if !(isnothing(H))
        # H[1,1] = ... ; H[2,2] = ... ; H[1,2] = ... ; H[2,1] = H[1,2]
    end
    if !(isnothing(F))
        #return cost
    end    
end

However, I get the following error message (starting from the optimize function above).

ERROR: MethodError: no method matching Cholesky{Float64, Matrix{Float64}}(::Cholesky{Float64, MMatrix{2, 2, Float64, 4}})
Closest candidates are:
  Cholesky{T, S}(::Any, ::Any, ::Any) where {T, S<:(AbstractMatrix)} at C:\Users\xx\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\LinearAlgebra\src\cholesky.jl:87
Stacktrace:
  [1] convert(#unused#::Type{Cholesky{Float64, Matrix{Float64}}}, f::Cholesky{Float64, MMatrix{2, 2, Float64, 4}})
    @ LinearAlgebra C:\Users\xx\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\LinearAlgebra\src\factorization.jl:56
  [2] setproperty!(x::Optim.NewtonState{MVector{2, Float64}, Float64, Cholesky{Float64, Matrix{Float64}}}, 
f::Symbol, v::Cholesky{Float64, MMatrix{2, 2, Float64, 4}})
    @ Base .\Base.jl:39
  [3] update_state!(d::TwiceDifferentiable{Float64, MVector{2, Float64}, MMatrix{2, 2, Float64, 4}, MVector{2, Float64}}, state::Optim.NewtonState{MVector{2, Float64}, Float64, Cholesky{Float64, Matrix{Float64}}}, 
method::Newton{LineSearches.InitialStatic{Float64}, LineSearches.HagerZhang{Float64, Base.RefValue{Bool}}})    @ Optim C:\Users\xx\.julia\packages\Optim\tP8PJ\src\multivariate\solvers\second_order\newton.jl:67  
  [4] optimize(d::TwiceDifferentiable{Float64, MVector{2, Float64}, MMatrix{2, 2, Float64, 4}, MVector{2, Float64}}, initial_x::MVector{2, Float64}, method::Newton{LineSearches.InitialStatic{Float64}, LineSearches.HagerZhang{Float64, Base.RefValue{Bool}}}, options::Optim.Options{Float64, Nothing}, state::Optim.NewtonState{MVector{2, Float64}, Float64, Cholesky{Float64, Matrix{Float64}}})
    @ Optim C:\Users\xx\.julia\packages\Optim\tP8PJ\src\multivariate\optimize\optimize.jl:54
  [5] optimize
    @ C:\Users\xx\.julia\packages\Optim\tP8PJ\src\multivariate\optimize\optimize.jl:36 [inlined]        
  [6] optimize(f::NLSolversBase.InplaceObjective{Nothing, Nothing, var"#176#177"{SVector{3, Float64}, SVector{3, Float64}, SVector{3, Float64}, SVector{3, Float64}, SVector{3, Float64}}, Nothing, Nothing}, initial_x::MVector{2, Float64}, method::Newton{LineSearches.InitialStatic{Float64}, LineSearches.HagerZhang{Float64, Base.RefValue{Bool}}}, options::Optim.Options{Float64, Nothing}; inplace::Bool, autodiff::Symbol)        
    @ Optim C:\Users\xx\.julia\packages\Optim\tP8PJ\src\multivariate\optimize\interface.jl:142
  [7] optimize(f::NLSolversBase.InplaceObjective{Nothing, Nothing, var"#176#177"{SVector{3, Float64}, SVector{3, Float64}, SVector{3, Float64}, SVector{3, Float64}, SVector{3, Float64}}, Nothing, Nothing}, initial_x::MVector{2, Float64}, method::Newton{LineSearches.InitialStatic{Float64}, LineSearches.HagerZhang{Float64, Base.RefValue{Bool}}}, options::Optim.Options{Float64, Nothing}) (repeats 2 times)
    @ Optim C:\Users\xx\.julia\packages\Optim\tP8PJ\src\multivariate\optimize\interface.jl:138

I also get the same error when I declare the types of F,G and H above (Float64,MVector,and MMatrix). In addition, a very similar error message appears when I use SArrays.

Thank you very much for your help!

[pkofod]

Yes, that is not supported. You can find SArray/MArray support in the attempt at rewriting Optim from scratch in https://github.com/JuliaNLSolvers/NLSolvers.jl that is also a registered package but less tested. Let me know if you need help.