Check primal feasibility after termination failure (#435)

src/MOI_wrapper.jl

@@ -1167,7 +1167,7 @@ function MOI.get(model::Optimizer, ::MOI.TerminationStatus)
     elseif status == :Maximum_WallTime_Exceeded
         return MOI.TIME_LIMIT
     elseif status == :Restoration_Failed
-        return MOI.NUMERICAL_ERROR
+        return MOI.OTHER_ERROR
     elseif status == :Error_In_Step_Computation
         return MOI.NUMERICAL_ERROR
     elseif status == :Invalid_Option
@@ -1202,6 +1202,31 @@ end
 
 ### MOI.PrimalStatus
 
+function _manually_evaluated_primal_status(model::Optimizer)
+    x, g = model.inner.x, model.inner.g
+    m, n = length(g), length(x)
+    x_L, x_U = model.variables.lower, model.variables.upper
+    g_L, g_U = copy(model.qp_data.g_L), copy(model.qp_data.g_U)
+    # Assuming constraints are guaranteed to be in the order [qp_cons, nlp_cons]
+    for bound in model.nlp_data.constraint_bounds
+        push!(g_L, bound.lower)
+        push!(g_U, bound.upper)
+    end
+    # 1e-8 is the default tolerance
+    tol = get(model.options, "tol", 1e-8)
+    if all(x_L[i] - tol <= x[i] <= x_U[i] + tol for i in 1:n) &&
+       all(g_L[i] - tol <= g[i] <= g_U[i] + tol for i in 1:m)
+        return MOI.FEASIBLE_POINT
+    end
+    # 1e-6 is the default acceptable tolerance
+    atol = get(model.options, "acceptable_tol", 1e-6)
+    if all(x_L[i] - atol <= x[i] <= x_U[i] + atol for i in 1:n) &&
+       all(g_L[i] - atol <= g[i] <= g_U[i] + atol for i in 1:m)
+        return MOI.NEARLY_FEASIBLE_POINT
+    end
+    return MOI.INFEASIBLE_POINT
+end
+
 function MOI.get(model::Optimizer, attr::MOI.PrimalStatus)
     if !(1 <= attr.result_index <= MOI.get(model, MOI.ResultCount()))
         return MOI.NO_SOLUTION
@@ -1218,7 +1243,9 @@ function MOI.get(model::Optimizer, attr::MOI.PrimalStatus)
     elseif status == :Infeasible_Problem_Detected
         return MOI.INFEASIBLE_POINT
     else
-        return MOI.UNKNOWN_RESULT_STATUS
+        # Not sure. RestorationFailure can terminate at a feasible (but
+        # non-stationary) point.
+        return _manually_evaluated_primal_status(model)
     end
 end
 

test/MOI_wrapper.jl

@@ -633,6 +633,45 @@ function test_nlp_model_set_set()
     return
 end
 
+function test_manually_evaluated_primal_status()
+    model = Ipopt.Optimizer()
+    MOI.set(model, MOI.Silent(), true)
+    x = MOI.add_variable(model)
+    MOI.add_constraint(model, x, MOI.Interval(0.0, 1.0))
+    MOI.add_constraint(model, 1.0 * x, MOI.Interval(0.1, 1.1))
+    f = MOI.ScalarNonlinearFunction(:log, Any[x])
+    MOI.add_constraint(model, f, MOI.Interval(-1.0, 1.0))
+    MOI.optimize!(model)
+    x_star = copy(model.inner.x)
+    g_star = copy(model.inner.g)
+    for (xi, status) in (
+        x_star[1] => MOI.FEASIBLE_POINT,
+        -1.0 => MOI.INFEASIBLE_POINT,
+        -1e-7 => MOI.NEARLY_FEASIBLE_POINT,
+        1 + 1e-7 => MOI.NEARLY_FEASIBLE_POINT,
+        0.0 => MOI.FEASIBLE_POINT,
+        1.0 => MOI.FEASIBLE_POINT,
+    )
+        model.inner.x[1] = xi
+        @test Ipopt._manually_evaluated_primal_status(model) == status
+    end
+    model.inner.x .= x_star
+    for (gi, status) in (
+        g_star => MOI.FEASIBLE_POINT,
+        [0.0, 0.0] => MOI.INFEASIBLE_POINT,
+        [0.1 - 1e-7, 0.0] => MOI.NEARLY_FEASIBLE_POINT,
+        [1.1 + 1e-7, 0.0] => MOI.NEARLY_FEASIBLE_POINT,
+        [0.1, -1.0 - 1e-7] => MOI.NEARLY_FEASIBLE_POINT,
+        [0.1, 1.0 + 1e-7] => MOI.NEARLY_FEASIBLE_POINT,
+        [0.1, 1.0] => MOI.FEASIBLE_POINT,
+        [1.1, -1.0] => MOI.FEASIBLE_POINT,
+    )
+        model.inner.g .= gi
+        @test Ipopt._manually_evaluated_primal_status(model) == status
+    end
+    return
+end
+
 end  # module TestMOIWrapper
 
 TestMOIWrapper.runtests()