lanl-ansi / Alpine.jl

A Julia/JuMP-based Global Optimization Solver for Non-convex Programs
https://lanl-ansi.github.io/Alpine.jl/latest/
Other
241 stars 39 forks source link

Alpine throws "There are currently 0 solution(s) in the model." for a solvable bilinear optimization problem #210

Closed Shuvomoy closed 1 year ago

Shuvomoy commented 1 year ago

I am trying to solve a bilinear optimization problem using Alpine. If I try Gurobi's bilinear solver it works fine, but Alpine does not.

Gurobi code

## Load packages
using Alpine
using JuMP
using Gurobi
using Ipopt
using KNITRO

## Solve using Gurobi

d = 10

nonlinear_model = Model(Gurobi.Optimizer)

set_optimizer_attribute(nonlinear_model, "MIPFocus", 3)

set_optimizer_attribute(nonlinear_model, "NonConvex", 2)

# declare the variable 
@variable(nonlinear_model, 0 <= x[1:d])

# add objective
@objective(nonlinear_model, Max, x[1] * x[d] - x[2] * x[d-1])

# add the constraints

@constraint(nonlinear_model, sum(x[i] for i in 1:d) <= 10)

for i in 1:d-1
    @constraint(nonlinear_model, x[i] * x[i+1] <= 2)
end

@constraint(nonlinear_model, sum(x[i] * x[i+1] for i in 1:d-1) == 1)

# optimize the model

optimize!(nonlinear_model)

Gurobi output is:

Set parameter Username
Academic license - for non-commercial use only - expires 2023-08-02
Set parameter MIPFocus to value 3
Set parameter NonConvex to value 2
Set parameter NonConvex to value 2
Set parameter MIPFocus to value 3
Gurobi Optimizer version 9.5.2 build v9.5.2rc0 (mac64[rosetta2])
Thread count: 10 physical cores, 10 logical processors, using up to 10 threads
Optimize a model with 1 rows, 10 columns and 10 nonzeros
Model fingerprint: 0x7855f5df
Model has 2 quadratic objective terms
Model has 10 quadratic constraints
Coefficient statistics:
  Matrix range     [1e+00, 1e+00]
  QMatrix range    [1e+00, 1e+00]
  Objective range  [0e+00, 0e+00]
  QObjective range [2e+00, 2e+00]
  Bounds range     [0e+00, 0e+00]
  RHS range        [1e+01, 1e+01]
  QRHS range       [1e+00, 2e+00]

Continuous model is non-convex -- solving as a MIP

Presolve time: 0.00s
Presolved: 61 rows, 22 columns, 120 nonzeros
Presolved model has 20 bilinear constraint(s)
Variable types: 22 continuous, 0 integer (0 binary)
Root relaxation presolve removed 39 rows and 2 columns
Root relaxation presolved: 22 rows, 20 columns, 59 nonzeros

Root relaxation: objective 4.955556e+01, 19 iterations, 0.00 seconds (0.00 work units)

    Nodes    |    Current Node    |     Objective Bounds      |     Work
 Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time

     0     0   49.55556    0   11          -   49.55556      -     -    0s
H    0     0                      23.9791576   49.55556   107%     -    0s
     0     0   24.99322    0   10   23.97916   24.99322  4.23%     -    0s
     0     0   24.99231    0    8   23.97916   24.99231  4.23%     -    0s
     0     0   24.99029    0   10   23.97916   24.99029  4.22%     -    0s
     0     0   24.99026    0   11   23.97916   24.99026  4.22%     -    0s
     0     2   24.99026    0   11   23.97916   24.99026  4.22%     -    0s
H    3     6                      23.9946925   24.95311  3.99%   5.3    0s
*   44    10               8      23.9965699   24.15915  0.68%   5.8    0s
*   50    12               9      23.9998195   24.02209  0.09%   5.3    0s
*   74    15              11      23.9999121   24.01252  0.05%   3.9    0s
H   76    15                      23.9999978   24.01252  0.05%   3.8    0s

Cutting planes:
  RLT: 18

Explored 147 nodes (367 simplex iterations) in 0.03 seconds (0.01 work units)
Thread count was 10 (of 10 available processors)

Solution count 6: 24 23.9999 23.9998 ... 23.9792

Optimal solution found (tolerance 1.00e-04)
Best objective 2.399999783339e+01, best bound 2.400193154381e+01, gap 0.0081%

User-callback calls 511, time in user-callback 0.00 sec

While solving the same problem using Alpoine, I am getting There are currently 0 solution(s) in the model.

Alpine code:

## Set MIP solver
gurobi = JuMP.optimizer_with_attributes(Gurobi.Optimizer, 
                                        MOI.Silent() => true, 
                                        "Presolve" => 2)

# Set NLP local solver
knitro = JuMP.optimizer_with_attributes(KNITRO.Optimizer, 
MOI.Silent() => true)

ipopt = JuMP.optimizer_with_attributes(Ipopt.Optimizer, 
MOI.Silent() => true)

# Set Alpine as the global solver
const alpine = JuMP.optimizer_with_attributes(Alpine.Optimizer, 
                                              "nlp_solver"   => ipopt, # knitro,  
                                              "mip_solver"   => gurobi,
                                              "presolve_bt"  => true,
                                              "disc_ratio"   => 10,
                                              "apply_partitioning" => true,
                                              "log_level" => 100
                                              )

m = nlp(solver = alpine)
JuMP.optimize!(m)

The output is:

PROBLEM STATISTICS
  Objective sense = Max
  # Variables = 10
  # Bin-Int Variables = 0
  # Constraints = 11
  # NL Constraints = 10
  # Linear Constraints = 1
  # Detected convex constraints = 0
  # Detected nonlinear terms = 11
  # Variables involved in nonlinear terms = 10
  # Potential variables for partitioning = 10
SUB-SOLVERS USED BY ALPINE
  NLP local solver = Ipopt
  MIP solver = Gurobi
ALPINE CONFIGURATION
  Maximum iterations (upper-bounding MIPs) =  99
  Relative global optimality gap = 0.01%
  Discretization ratio = 10
  Bound-tightening presolve = true
  Maximum iterations (OBBT) = 25
PRESOLVE 
  Doing local search

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
 Ipopt is released as open source code under the Eclipse Public License (EPL).
         For more information visit https://github.com/coin-or/Ipopt
******************************************************************************

  Local solver returns a feasible point with value 23.9792
  Starting bound-tightening
Set parameter Username
Academic license - for non-commercial use only - expires 2023-08-02
ERROR: Result index of attribute MathOptInterface.ObjectiveValue(1) out of bounds. There are currently 0 solution(s) in the model.
Stacktrace:
  [1] check_result_index_bounds
    @ ~/.julia/packages/MathOptInterface/LqT2k/src/attributes.jl:198 [inlined]
  [2] get(model::Gurobi.Optimizer, attr::MathOptInterface.ObjectiveValue)
    @ Gurobi ~/.julia/packages/Gurobi/VPomg/src/MOI_wrapper/MOI_wrapper.jl:3084
  [3] get(b::MathOptInterface.Bridges.LazyBridgeOptimizer{Gurobi.Optimizer}, attr::MathOptInterface.ObjectiveValue)
    @ MathOptInterface.Bridges ~/.julia/packages/MathOptInterface/LqT2k/src/Bridges/bridge_optimizer.jl:988
  [4] _get_model_attribute(model::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.Bridges.LazyBridgeOptimizer{Gurobi.Optimizer}, MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}}, attr::MathOptInterface.ObjectiveValue)
    @ MathOptInterface.Utilities ~/.julia/packages/MathOptInterface/LqT2k/src/Utilities/cachingoptimizer.jl:828
  [5] get
    @ ~/.julia/packages/MathOptInterface/LqT2k/src/Utilities/cachingoptimizer.jl:876 [inlined]
  [6] _moi_get_result(model::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.Bridges.LazyBridgeOptimizer{Gurobi.Optimizer}, MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}}, args::MathOptInterface.ObjectiveValue)
    @ JuMP ~/.julia/packages/JuMP/60Bnj/src/JuMP.jl:1155
  [7] get(model::Model, attr::MathOptInterface.ObjectiveValue)
    @ JuMP ~/.julia/packages/JuMP/60Bnj/src/JuMP.jl:1175
  [8] objective_value(model::Model; result::Int64)
    @ JuMP ~/.julia/packages/JuMP/60Bnj/src/objective.jl:42
  [9] objective_value
    @ ~/.julia/packages/JuMP/60Bnj/src/objective.jl:41 [inlined]
 [10] solve_bound_tightening_model(m::Alpine.Optimizer; kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ Alpine ~/.julia/packages/Alpine/RIjEU/src/presolve.jl:285
 [11] solve_bound_tightening_model
    @ ~/.julia/packages/Alpine/RIjEU/src/presolve.jl:265 [inlined]
 [12] minmax_bound_tightening(m::Alpine.Optimizer; use_bound::Bool, time_limit::Float64, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ Alpine ~/.julia/packages/Alpine/RIjEU/src/presolve.jl:110
 [13] bound_tightening(m::Alpine.Optimizer; use_bound::Bool, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ Alpine ~/.julia/packages/Alpine/RIjEU/src/presolve.jl:18
 [14] presolve(m::Alpine.Optimizer)
    @ Alpine ~/.julia/packages/Alpine/RIjEU/src/algorithm.jl:201
 [15] optimize!(m::Alpine.Optimizer)
    @ Alpine ~/.julia/packages/Alpine/RIjEU/src/algorithm.jl:143
 [16] optimize!
    @ ~/.julia/packages/MathOptInterface/LqT2k/src/Bridges/bridge_optimizer.jl:376 [inlined]
 [17] optimize!
    @ ~/.julia/packages/MathOptInterface/LqT2k/src/MathOptInterface.jl:87 [inlined]
 [18] optimize!(m::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.Bridges.LazyBridgeOptimizer{Alpine.Optimizer}, MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.Model{Float64}}})
    @ MathOptInterface.Utilities ~/.julia/packages/MathOptInterface/LqT2k/src/Utilities/cachingoptimizer.jl:316
 [19] optimize!(model::Model; ignore_optimize_hook::Bool, _differentiation_backend::MathOptInterface.Nonlinear.SparseReverseMode, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ JuMP ~/.julia/packages/JuMP/60Bnj/src/optimizer_interface.jl:185
 [20] optimize!(model::Model)
    @ JuMP ~/.julia/packages/JuMP/60Bnj/src/optimizer_interface.jl:155
 [21] top-level scope
    @ ~/Desktop/alpine_test.jl:101

(updated the issue as suggested by @blegat )

harshangrjn commented 1 year ago

@Shuvomoy @blegat I suspect the issue is not with Alpine, and rather is in Gurobi's presolve while solving linear programs (in your case) within Alpine's OBBT presolve. This seems to be an issue which we have observed on many other QCQPs. Here is the output below with "Presolve" => 0 while invoking the Gurobi solver and other settings of Alpine as you showed above. Also tested with CPLEX, with both it's presolve turned on and off. In both settings, Alpine+CPLEX seems to converge without any issue.

PROBLEM STATISTICS
  Objective sense = Max
  # Variables = 10
  # Bin-Int Variables = 0
  # Constraints = 11
  # NL Constraints = 10
  # Linear Constraints = 1
  # Detected convex constraints = 0
  # Detected nonlinear terms = 11
  # Variables involved in nonlinear terms = 10
  # Potential variables for partitioning = 10
SUB-SOLVERS USED BY ALPINE
  NLP local solver = Ipopt
  MIP solver = Gurobi
ALPINE CONFIGURATION
  Maximum iterations (upper-bounding MIPs) =  99
  Relative global optimality gap = 0.01%
  Discretization ratio = 10
  Bound-tightening presolve = true
  Maximum iterations (OBBT) = 25
PRESOLVE 
  Doing local search
  Local solver returns a feasible point with value 23.9792
  Starting bound-tightening
  Actual iterations (OBBT): 9
  Post-presolve optimality gap: 0.569%
  Completed presolve in 0.57s
UPPER-BOUNDING ITERATIONS
====================================================================================================
| Iter   | Incumbent       | Best Incumbent      | Upper Bound        | Gap (%)         | Time      
| 1      | 23.9792         | 23.9792             | 24.009             | 0.124           | 0.59s            
| 2      | 24.0            | 24.0                | 24.009             | 0.038           | 0.65s            
| finish | 24.0            | 24.0                | 24.0021            | 0.009           | 0.77s            
====================================================================================================

*** Alpine ended with status OPTIMAL ***

Of course, Alpine can also be run without it's OBBT presolve turned on (setting as "presolve_bt" => false). By doing so, it converges fine with Gurobi's presolve on. Do let me know if in case you observe more issues.

Shuvomoy commented 1 year ago

Okay I see thanks @harshangrjn !