diff --git a/ext/MadNLPMOI/MadNLPMOI.jl b/ext/MadNLPMOI/MadNLPMOI.jl index f9c7234b..7467e16b 100644 --- a/ext/MadNLPMOI/MadNLPMOI.jl +++ b/ext/MadNLPMOI/MadNLPMOI.jl @@ -99,6 +99,7 @@ MOI.eval_constraint(::_EmptyNLPEvaluator, g, x) = nothing MOI.jacobian_structure(::_EmptyNLPEvaluator) = Tuple{Int64,Int64}[] MOI.hessian_lagrangian_structure(::_EmptyNLPEvaluator) = Tuple{Int64,Int64}[] MOI.eval_constraint_jacobian(::_EmptyNLPEvaluator, J, x) = nothing +MOI.eval_constraint_jacobian_transpose_product(::_EmptyNLPEvaluator, Jtv, x, v) = nothing MOI.eval_hessian_lagrangian(::_EmptyNLPEvaluator, H, x, σ, μ) = nothing function MOI.empty!(model::Optimizer) @@ -123,6 +124,9 @@ function MOI.empty!(model::Optimizer) if haskey(model.options, :jacobian_constant) delete!(model.options, :jacobian_constant) end + if haskey(model.options, :hessian_approximation) + delete!(model.options, :hessian_approximation) + end return end @@ -141,6 +145,16 @@ function MOI.copy_to(model::Optimizer, src::MOI.ModelLike) return MOI.Utilities.default_copy_to(model, src) end +function _init_nlp_model(model) + if model.nlp_model === nothing + if !(model.nlp_data.evaluator isa _EmptyNLPEvaluator) + error("Cannot mix the new and legacy nonlinear APIs") + end + model.nlp_model = MOI.Nonlinear.Model() + end + return +end + MOI.get(::Optimizer, ::MOI.SolverName) = "MadNLP" function MOI.supports_add_constrained_variable( @@ -487,6 +501,28 @@ function MOI.get(model::Optimizer, attr::MOI.ListOfSupportedNonlinearOperators) return MOI.get(model.nlp_model, attr) end +### UserDefinedFunction + +MOI.supports(model::Optimizer, ::MOI.UserDefinedFunction) = true + +function MOI.set(model::Optimizer, attr::MOI.UserDefinedFunction, args) + _init_nlp_model(model) + MOI.Nonlinear.register_operator( + model.nlp_model, + attr.name, + attr.arity, + args..., + ) + return +end + +### ListOfSupportedNonlinearOperators + +function MOI.get(model::Optimizer, attr::MOI.ListOfSupportedNonlinearOperators) + _init_nlp_model(model) + return MOI.get(model.nlp_model, attr) +end + ### MOI.VariablePrimalStart function MOI.supports( @@ -726,6 +762,13 @@ function MOI.eval_constraint_jacobian(model::Optimizer, values, x) return end +function MOI.eval_constraint_jacobian_transpose_product(model::Optimizer, Jtv, x, v) + MOI.eval_constraint_jacobian_transpose_product(model.nlp_data.evaluator, Jtv, x, v) + # Evaluate QPBlockData after NLPEvaluator to ensure that Jtv is not reset. + MOI.eval_constraint_jacobian_transpose_product(model.qp_data, Jtv, x, v) + return +end + ### Eval_H_CB function MOI.hessian_lagrangian_structure(model::Optimizer) @@ -751,23 +794,27 @@ end NLPModels.obj(nlp::MOIModel, x::AbstractVector{Float64}) = MOI.eval_objective(nlp.model,x) -function NLPModels. grad!(nlp::MOIModel, x::AbstractVector{Float64}, g::AbstractVector{Float64}) +function NLPModels.grad!(nlp::MOIModel, x::AbstractVector{Float64}, g::AbstractVector{Float64}) MOI.eval_objective_gradient(nlp.model, g, x) end -function NLPModels. cons!(nlp::MOIModel, x::AbstractVector{Float64}, c::AbstractVector{Float64}) +function NLPModels.cons!(nlp::MOIModel, x::AbstractVector{Float64}, c::AbstractVector{Float64}) MOI.eval_constraint(nlp.model, c, x) end -function NLPModels. jac_coord!(nlp::MOIModel, x::AbstractVector{Float64}, jac::AbstractVector{Float64}) +function NLPModels.jac_coord!(nlp::MOIModel, x::AbstractVector{Float64}, jac::AbstractVector{Float64}) MOI.eval_constraint_jacobian(nlp.model, jac, x) end -function NLPModels. hess_coord!(nlp::MOIModel, x::AbstractVector{Float64}, l::AbstractVector{Float64}, hess::AbstractVector{Float64}; obj_weight::Float64=1.0) +function NLPModels.jtprod!(nlp::MOIModel, x::AbstractVector{Float64}, v::Vector{Float64}, Jtv::AbstractVector{Float64}) + MOI.eval_constraint_jacobian_transpose_product(nlp.model, Jtv, x, v) +end + +function NLPModels.hess_coord!(nlp::MOIModel, x::AbstractVector{Float64}, l::AbstractVector{Float64}, hess::AbstractVector{Float64}; obj_weight::Float64=1.0) MOI.eval_hessian_lagrangian(nlp.model, hess, x, obj_weight, l) end -function NLPModels. hess_structure!(nlp::MOIModel, I::AbstractVector{T}, J::AbstractVector{T}) where T +function NLPModels.hess_structure!(nlp::MOIModel, I::AbstractVector{T}, J::AbstractVector{T}) where T @assert length(I) == length(J) == length(MOI.hessian_lagrangian_structure(nlp.model)) cnt = 1 for (row, col) in MOI.hessian_lagrangian_structure(nlp.model) @@ -804,11 +851,11 @@ function MOIModel(model::Optimizer) any(isequal(_kFunctionTypeScalarQuadratic), model.qp_data.function_type) has_nlp_constraints = !isempty(model.nlp_data.constraint_bounds) has_nlp_objective = model.nlp_data.has_objective - has_nlp_hessian = :Hess in MOI.features_available(model.nlp_data.evaluator) + has_hessian = :Hess in MOI.features_available(model.nlp_data.evaluator) is_nlp = has_nlp_constraints || has_nlp_objective # Initialize evaluator using model's structure. init_feat = [:Grad] - if has_nlp_hessian + if has_hessian push!(init_feat, :Hess) end if has_nlp_constraints @@ -836,7 +883,11 @@ function MOIModel(model::Optimizer) # Sparsity jacobian_sparsity = MOI.jacobian_structure(model) - hessian_sparsity = MOI.hessian_lagrangian_structure(model) + hessian_sparsity = if has_hessian + MOI.hessian_lagrangian_structure(model) + else + Tuple{Int,Int}[] + end nnzh = length(hessian_sparsity) nnzj = length(jacobian_sparsity) @@ -858,6 +909,9 @@ function MOIModel(model::Optimizer) if !has_nlp_constraints && !has_quadratic_constraints model.options[:jacobian_constant] = true end + if !has_hessian + model.options[:hessian_approximation] = MadNLP.CompactLBFGS + end return MOIModel( NLPModels.NLPModelMeta( @@ -954,7 +1008,7 @@ function MOI.get(model::Optimizer, ::MOI.RawStatusString) elseif model.solver === nothing return "Optimize not called" end - return MadNLP.get_status_output(model.result.status, model.result.options) + return MadNLP.get_status_output(model.result.status, model.result.options) end diff --git a/ext/MadNLPMOI/utils.jl b/ext/MadNLPMOI/utils.jl index 9c703719..17433635 100644 --- a/ext/MadNLPMOI/utils.jl +++ b/ext/MadNLPMOI/utils.jl @@ -201,24 +201,25 @@ function eval_sparse_gradient( f::MOI.ScalarQuadraticFunction{T}, x::Vector{T}, p::Dict{Int64,T}, + adj::T, )::Int where {T} i = 0 for term in f.affine_terms if !_is_parameter(term.variable) i += 1 - ∇f[i] = term.coefficient + ∇f[i] += term.coefficient * adj end end for term in f.quadratic_terms if !_is_parameter(term.variable_1) v = _value(term.variable_2, x, p) i += 1 - ∇f[i] = term.coefficient * v + ∇f[i] += term.coefficient * v * adj end if term.variable_1 != term.variable_2 && !_is_parameter(term.variable_2) v = _value(term.variable_1, x, p) i += 1 - ∇f[i] = term.coefficient * v + ∇f[i] += term.coefficient * v * adj end end return i @@ -229,12 +230,13 @@ function eval_sparse_gradient( f::MOI.ScalarAffineFunction{T}, x::Vector{T}, p::Dict{Int64,T}, + adj::T, )::Int where {T} i = 0 for term in f.terms if !_is_parameter(term.variable) i += 1 - ∇f[i] = term.coefficient + ∇f[i] += term.coefficient * adj end end return i @@ -508,13 +510,26 @@ function MOI.eval_constraint_jacobian( x::AbstractVector{T}, ) where {T} i = 1 + fill!(J, zero(T)) for constraint in block.constraints ∇f = view(J, i:length(J)) - i += eval_sparse_gradient(∇f, constraint, x, block.parameters) + i += eval_sparse_gradient(∇f, constraint, x, block.parameters, one(T)) end return i end +function MOI.eval_constraint_jacobian_transpose_product( + block::QPBlockData{T}, + Jtv::AbstractVector{T}, + x::AbstractVector{T}, + v::AbstractVector{T}, +) where {T} + for (i, constraint) in enumerate(block.constraints) + eval_sparse_gradient(∇f, constraint, x, block.parameters, v[i]) + end + return Jtv +end + function MOI.hessian_lagrangian_structure(block::QPBlockData) H = Tuple{Int,Int}[] append_sparse_hessian_structure!(block.objective, H) diff --git a/lib/MadNLPTests/src/Instances/hs15nohessian.jl b/lib/MadNLPTests/src/Instances/hs15nohessian.jl new file mode 100644 index 00000000..52a2021f --- /dev/null +++ b/lib/MadNLPTests/src/Instances/hs15nohessian.jl @@ -0,0 +1,65 @@ +struct HS15NoHessianModel{T} <: NLPModels.AbstractNLPModel{T,Vector{T}} + meta::NLPModels.NLPModelMeta{T, Vector{T}} + counters::NLPModels.Counters +end + +function HS15NoHessianModel(;T = Float64, x0=zeros(T,2), y0=zeros(T,2)) + return HS15NoHessianModel( + NLPModels.NLPModelMeta( + 2, #nvar + ncon = 2, + nnzj = 4, + nnzh = 0, + x0 = x0, + y0 = y0, + lvar = T[-Inf, -Inf], + uvar = T[0.5, Inf], + lcon = T[1.0, 0.0], + ucon = T[Inf, Inf], + minimize = true + ), + NLPModels.Counters() + ) +end + +function NLPModels.obj(nlp::HS15NoHessianModel, x::AbstractVector) + return 100.0 * (x[2] - x[1]^2)^2 + (1.0 - x[1])^2 +end + +function NLPModels.grad!(nlp::HS15NoHessianModel, x::AbstractVector, g::AbstractVector) + z = x[2] - x[1]^2 + g[1] = -400.0 * z * x[1] - 2.0 * (1.0 - x[1]) + g[2] = 200.0 * z + return +end + +function NLPModels.cons!(nlp::HS15NoHessianModel, x::AbstractVector, c::AbstractVector) + c[1] = x[1] * x[2] + c[2] = x[1] + x[2]^2 +end + +function NLPModels.jac_structure!(nlp::HS15NoHessianModel, I::AbstractVector{T}, J::AbstractVector{T}) where T + copyto!(I, [1, 1, 2, 2]) + copyto!(J, [1, 2, 1, 2]) +end + +function NLPModels.jac_coord!(nlp::HS15NoHessianModel, x::AbstractVector, J::AbstractVector) + J[1] = x[2] # (1, 1) + J[2] = x[1] # (1, 2) + J[3] = 1.0 # (2, 1) + J[4] = 2*x[2] # (2, 2) + return J +end + +function NLPModels.jprod!(nlp::HS15NoHessianModel, x::AbstractVector, v::AbstractVector, jv::AbstractVector) + jv[1] = x[2] * v[1] + x[1] * v[2] + jv[2] = v[1] + 2 * x[2] * v[2] + return jv +end + +function NLPModels.jtprod!(nlp::HS15NoHessianModel, x::AbstractVector, v::AbstractVector, jv::AbstractVector) + jv[1] = x[2] * v[1] + v[2] + jv[2] = x[1] * v[1] + 2 * x[2] * v[2] + return jv +end + diff --git a/lib/MadNLPTests/src/MadNLPTests.jl b/lib/MadNLPTests/src/MadNLPTests.jl index 7bf26e43..b5012fdf 100644 --- a/lib/MadNLPTests/src/MadNLPTests.jl +++ b/lib/MadNLPTests/src/MadNLPTests.jl @@ -354,6 +354,7 @@ end include("Instances/dummy_qp.jl") include("Instances/hs15.jl") +include("Instances/hs15nohessian.jl") include("Instances/nls.jl") include("wrapper.jl") diff --git a/src/nlpmodels.jl b/src/nlpmodels.jl index b4b773b9..8336995a 100644 --- a/src/nlpmodels.jl +++ b/src/nlpmodels.jl @@ -357,7 +357,9 @@ function create_callback( jac_scale = similar(jac_buffer, nnzj) ; fill!(jac_scale, one(T)) NLPModels.jac_structure!(nlp, jac_I, jac_J) - NLPModels.hess_structure!(nlp, hess_I, hess_J) + if nnzh > 0 + NLPModels.hess_structure!(nlp, hess_I, hess_J) + end fixed_handler, nnzj, nnzh = create_sparse_fixed_handler( fixed_variable_treatment, diff --git a/test/MOI_interface_test.jl b/test/MOI_interface_test.jl index 6151d96d..750e4fc1 100644 --- a/test/MOI_interface_test.jl +++ b/test/MOI_interface_test.jl @@ -50,8 +50,6 @@ function test_MOI_Test() # - Excluded because Hessian information is needed "test_nonlinear_hs071_hessian_vector_product", # - Excluded because Hessian information is needed - "test_nonlinear_hs071_no_hessian", - # - Excluded because Hessian information is needed "test_nonlinear_invalid", # - Excluded because this test is optional @@ -99,6 +97,20 @@ function test_invalid_number_in_hessian_lagrangian() return end +# See issue #318 +function test_user_defined_function() + model = MadNLP.Optimizer() + MOI.set(model, MOI.Silent(), true) + # Define custom function. + f(a, b) = a^2 + b^2 + x = MOI.add_variables(model, 2) + MOI.set(model, MOI.UserDefinedFunction(:f, 2), (f,)) + obj_f = MOI.ScalarNonlinearFunction(:f, Any[x[1], x[2]]) + MOI.set(model, MOI.ObjectiveFunction{typeof(obj_f)}(), obj_f) + MOI.optimize!(model) + @test MOI.get(model, MOI.TerminationStatus()) == MOI.LOCALLY_SOLVED +end + end TestMOIWrapper.runtests() diff --git a/test/madnlp_quasi_newton.jl b/test/madnlp_quasi_newton.jl index 1f916d93..88cbb001 100644 --- a/test/madnlp_quasi_newton.jl +++ b/test/madnlp_quasi_newton.jl @@ -34,7 +34,7 @@ end @testset "MadNLP: LBFGS" begin @testset "HS15" begin - nlp = MadNLPTests.HS15Model() + nlp = MadNLPTests.HS15NoHessianModel() solver_qn = MadNLP.MadNLPSolver( nlp; callback = MadNLP.SparseCallback,