search

jump-dev / HiGHS.jl

A Julia interface to the HiGHS solver
https://highs.dev
MIT License
112 stars 16 forks source link

Getting error when retrieving results of optimization #60

Closed SupplyChef closed 3 years ago

SupplyChef commented 3 years ago

model.zip

With the model attached to this, if I run 

julia> m = read_from_file("model.mps")
julia> set_optimizer(m, Cbc.Optimizer)
julia> optimize!(m)
julia> termination_status(m)
julia> objective_value(m)

I get:

OPTIMAL::TerminationStatusCode = 1
5.234147524505995e8

When I do the same with HiGHS (i.e, set_optimizer(m, HiGHS.Optimizer) I get:

OPTIMAL::TerminationStatusCode = 1

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
   @ C:\Users\renau\.julia\packages\MathOptInterface\YDdD3\src\attributes.jl:139 [inlined]
 [2] get
   @ C:\Users\renau\.julia\packages\HiGHS\UdLoV\src\MOI_wrapper.jl:1669 [inlined]
 [3] get(b::MathOptInterface.Bridges.LazyBridgeOptimizer{HiGHS.Optimizer}, attr::MathOptInterface.ObjectiveValue)
   @ MathOptInterface.Bridges C:\Users\renau\.julia\packages\MathOptInterface\YDdD3\src\Bridges\bridge_optimizer.jl:913
 [4] get(model::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.AbstractOptimizer, MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.GenericModel{Float64, MathOptInterface.Utilities.ModelFunctionConstraints{Float64}}}}, attr::MathOptInterface.ObjectiveValue)
   @ MathOptInterface.Utilities C:\Users\renau\.julia\packages\MathOptInterface\YDdD3\src\Utilities\cachingoptimizer.jl:716
 [5] _moi_get_result(model::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.AbstractOptimizer, MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.GenericModel{Float64, MathOptInterface.Utilities.ModelFunctionConstraints{Float64}}}}, args::MathOptInterface.ObjectiveValue)
   @ JuMP C:\Users\renau\.julia\packages\JuMP\b3hGi\src\JuMP.jl:1199
 [6] get(model::Model, attr::MathOptInterface.ObjectiveValue)
   @ JuMP C:\Users\renau\.julia\packages\JuMP\b3hGi\src\JuMP.jl:1212
 [7] objective_value(model::Model; result::Int64)
   @ JuMP C:\Users\renau\.julia\packages\JuMP\b3hGi\src\objective.jl:42
 [8] objective_value(model::Model)
   @ JuMP C:\Users\renau\.julia\packages\JuMP\b3hGi\src\objective.jl:42
 [9] top-level scope
   @ REPL[32]:1

The screen output shows that HiGHS found the optimal solution.

odow commented 3 years ago

For some reason

julia> primal_status(m)
NO_SOLUTION::ResultStatusCode = 0

I'm guessing HiGHS changed: https://github.com/jump-dev/HiGHS.jl/blob/94427d21b67d55efe1dc6187144777379fd612fe/src/MOI_wrapper.jl#L1586-L1587

Here's what I get

julia> Highs_getIntInfoValue(model, "primal_solution_status", statusP)
0

julia> statusP
Base.RefValue{Int32}(1)

julia> Highs_getIntInfoValue(model, "dual_solution_status", statusP)
0

julia> statusP
Base.RefValue{Int32}(0)

Which means that HiGHS converged to an infeasible point: https://github.com/ERGO-Code/HiGHS/blob/d6a134805cfd7a057752c4f8b683bdeb7adf4bda/src/lp_data/HighsInfo.h#L224-L228

You can also see this in the report:

Solving report
  Status            Optimal
  Primal bound      523106468.392
  Dual bound        523106468.392
  Solution status   feasible
                    523106468.392 (objective)
                    9.09494701773e-13 (bound viol.)
                    0 (int. viol.)
                    1.87079422176e-07 (row viol.)
  Timing            3.24 (total)
                    0.35 (presolve)
                    0.00 (postsolve)
  Nodes             1
  LP iterations     867 (total)
                    0 (strong br.)
                    67 (separation)
                    0 (heuristics)

The row violation is 1.8e-07, which I guess exceeds the tolerances. You could try resolving this with reduced tolerances.

julia> set_optimizer_attribute(m, "primal_feasibility_tolerance", 1e-5)
1.0e-5

julia> optimize!(m)

Presolving model
950 rows, 37600 cols, 150200 nonzeros
Sparsify removed 0.0% of nonzeros
850 rows, 37600 cols, 112600 nonzeros
50 probing evaluations: 0 deleted rows, 0 deleted columns, 0 lifted nonzeros
850 rows, 37600 cols, 75150 nonzeros

Solving root node LP relaxation

Presolving model
850 rows, 37600 cols, 75150 nonzeros
800 rows, 37550 cols, 75050 nonzeros
800 rows, 37550 cols, 75050 nonzeros
Presolve : Reductions: rows 800(-50); columns 37550(-50); elements 75050(-100)
Solving the presolved LP
WARNING: Number of OMP threads available = 0 < 1 = Number of HiGHS threads to be used: Parallel performance will be less than anticipated
Using EKK dual simplex solver - serial
       Iteration        Objective     Infeasibilities num(sum)
DuPh2          0     0.0000000000e+00 Pr: 750(4458.07)
DuPh2         78     9.6783733134e+07 Pr: 720(2781.97)
DuPh2        155     1.7301577925e+08 Pr: 644(2057.71)
DuPh2        244     2.5680936447e+08 Pr: 555(1527.67)
DuPh2        339     3.1417623884e+08 Pr: 461(1140.94)
DuPh2        440     3.6373454461e+08 Pr: 360(811.272)
DuPh2        547     4.1385446769e+08 Pr: 253(524.21)
DuPh2        660     4.5005631050e+08 Pr: 140(267.826)
DuPh2        779     4.8748912711e+08 Pr: 21(36.8939)
DuPh2        800     4.9482727555e+08 Pr: 0(0)
DuPh2        800     4.9480910576e+08 Pr: 0(0)
Solving the original LP from the solution after postsolve
      time | open nodes |      nodes |     leaves |    lpiters | dual bound     | primal bound   | cutpool |  lpcuts |      gap | explored
 S    0.5s |          0 |          0 |          0 |          0 | -inf           | 523791824      |       0 |       0 |     inf% |    0.00%
 R    0.6s |          0 |          0 |          0 |        800 | 494809106      | 523791658      |       0 |       0 |    5.53% |    0.00%
      0.6s |          0 |          0 |          0 |        800 | 494809106      | 523791658      |       0 |       0 |    5.53% |    0.00%
      1.2s |          0 |          0 |          0 |        848 | 519742878      | 523791658      |     750 |      48 |    0.77% |    0.00%
      1.6s |          0 |          0 |          0 |        855 | 520655346      | 523791658      |     755 |      55 |    0.60% |    0.00%
      1.8s |          0 |          0 |          0 |        859 | 521132344      | 523791658      |     760 |      59 |    0.51% |    0.00%
      2.1s |          0 |          0 |          0 |        861 | 521220257      | 523791658      |     761 |      61 |    0.49% |    0.00%
 S    2.4s |          0 |          0 |          0 |        861 | 521220257      | 523578783      |     763 |      63 |    0.45% |    0.00%
      2.7s |          0 |          0 |          0 |        863 | 521378770      | 523578783      |     763 |      63 |    0.42% |    0.00%
      2.9s |          0 |          0 |          0 |        865 | 521909420      | 523578783      |     763 |      65 |    0.32% |    0.00%
 T    3.0s |          0 |          0 |          0 |        867 | 523106468      | 523106468      |     763 |      67 |    0.00% |    0.00%

Solving report
  Status            Optimal
  Primal bound      523106468.392
  Dual bound        523106468.392
  Solution status   feasible
                    523106468.392 (objective)
                    9.09494701773e-13 (bound viol.)
                    0 (int. viol.)
                    1.87079422176e-07 (row viol.)
  Timing            3.00 (total)
                    0.34 (presolve)
                    0.00 (postsolve)
  Nodes             1
  LP iterations     867 (total)
                    0 (strong br.)
                    67 (separation)
                    0 (heuristics)

julia> primal_status(m)
FEASIBLE_POINT::ResultStatusCode = 1

julia> objective_value(m)
5.231064683921523e8

This is not a bug in HiGHS.jl. You could open an issue at https://github.com/ERGO-Code/HiGHS/issues (with the same MPS file and text above) and link back to this one. They should probably say something other than "Optimal" and "feasible" in this case, since the solution clearly isn't.

SupplyChef commented 3 years ago

Thank you for the analysis. I opened an issue for HiGHS (https://github.com/ERGO-Code/HiGHS/issues/559).

lgottwald commented 3 years ago

The solution is feasible if the default feasibility tolerance for MIP of 1e-6 was not changed. The "primal_feasibility_tolerance" options is for LP and the "mip_feasibility_tolerance" option is for MIP.

lgottwald commented 3 years ago

I see, however, that there may be a bug in HiGHS that uses the LP feasibility tolerance for the MIP solution. I will check if this is the case.

odow commented 3 years ago

Closing since this is upstream in HiGHS