search

ERGO-Code / HiGHS

Linear optimization software
MIT License
935 stars 174 forks source link

HiGHS does not seem to terminate on AUG2DC problem #992

Open stephane-caron opened 1 year ago

stephane-caron commented 1 year ago

While testing HiGHS 1.2.2 (date: 2022-09-04, git hash: 8701dbf19) on the Maros-Meszaros test set, I found that HiGHS does not terminate on the AUG2DC problem (despite setting the time limit to 10 seconds).

Reproduction steps

Decompress AUG2DC.mat from AUG2DC.zip, then run the following Python code:

import highspy
import numpy as np
import scipy.io as spio
import scipy.sparse as spa

m = spio.loadmat("AUG2DC.mat", squeeze_me=True)
P = m["P"].astype(float).tocsc()
q = m["q"].astype(float)
G = m["G"].astype(float).tocsc()
h = m["h"].astype(float)
A = m["A"].astype(float).tocsc()
b = m["b"].astype(float)
lb = m["lb"].astype(float)
ub = m["ub"].astype(float)
model = highspy.HighsModel()

# Hessian
model.hessian_.dim_ = P.shape[0]
model.hessian_.start_ = P.indptr
model.hessian_.index_ = P.indices
model.hessian_.value_ = P.data

# Columns
n = q.shape[0]
model.lp_.num_col_ = n
model.lp_.col_cost_ = q
model.lp_.col_lower_ = lb
model.lp_.col_upper_ = ub

# Rows
model.lp_.num_row_ = 0
row_list: list = []
row_lower: list = []
row_upper: list = []
model.lp_.num_row_ += G.shape[0]
row_list.append(G)
row_lower.append(np.full((G.shape[0],), -highspy.kHighsInf))
row_upper.append(h)
model.lp_.num_row_ += A.shape[0]
row_list.append(A)
row_lower.append(b)
row_upper.append(b)
row_matrix = spa.vstack(row_list, format="csc")
model.lp_.a_matrix_.format_ = highspy.MatrixFormat.kColwise
model.lp_.a_matrix_.start_ = row_matrix.indptr
model.lp_.a_matrix_.index_ = row_matrix.indices
model.lp_.a_matrix_.value_ = row_matrix.data
model.lp_.a_matrix_.num_row_ = row_matrix.shape[0]
model.lp_.a_matrix_.num_col_ = row_matrix.shape[1]
model.lp_.row_lower_ = np.hstack(row_lower)
model.lp_.row_upper_ = np.hstack(row_upper)

# Call solver
solver = highspy.Highs()
solver.setOptionValue("log_to_console", True)
solver.setOptionValue("log_dev_level", highspy.HighsLogType.kVerbose)
solver.setOptionValue("highs_debug_level", highspy.HighsLogType.kVerbose)
solver.setOptionValue("time_limit", 10.0)
print("The solver should take no more than 10.0 seconds...")
solver.passModel(model)
solver.run()

The output on my machine is:

The solver should take no more than 10.0 seconds...
Running HiGHS 1.2.2 [date: 2022-09-04, git hash: 8701dbf19]
Copyright (c) 2022 ERGO-Code under MIT licence terms
Cols:       20200 lower bounds exceeding       -1e+20 are treated as -Infinity
Cols:       20200 upper bounds exceeding        1e+20 are treated as +Infinity
Hessian has dimension 20200 and 20200 nonzeros: inserting 0 zeros onto the diagonal
Running with 2 thread(s)
checkOptions: Options are OK
0, 2158350.216835, 10200, 0.564013, 0.000000, 0, 0.000000, 0.000000

Then the solver just hangs for at least 30 min.

Other details

I installed highspy==1.1.2.dev3 from PyPI.

feldmeier commented 1 year ago

Thanks for reporting this issues. The problem is that there is a 10200x10200 matrix that needs to be factorized right at the start of the first iteration, and the time limit is not checked while doing this (only at the end of every iteration). Chances are this factorization eats too much memory and forces swapping (and hence won't finish any time soon).

Dealing with situations like this (few active constraints at start point -> large initial nullspace -> large dense matrix to factorize) is a known bottleneck for the Nullspace Active Set QP Solver in HiGHS. I'm open to any suggestions to improve this behaviour, but chances are its inherent to the algorithm implemented here.

I guess it wouldn't hurt to add occasional checks for the time limit in some more parts of the qp solver code.

I would be happy to hear more about your experience using Highs on maros-meszaros! :-)

stephane-caron commented 1 year ago

Sure! I'm running it right now :slightly_smiling_face: The benchmark code is at https://github.com/stephane-caron/qpsolvers_benchmark and to running HiGHS on it is done by:

$ python run_benchmark.py run -s highs [-p PROBLEM_ID]

So far the solver has complained on the following problems:

stephane-caron commented 1 year ago

I would be happy to hear more about your experience using Highs on maros-meszaros! :-)

Results with HiGHS 1.2.2 and default solver settings:

More details in this report and the corresponding CSV data.

This is ongoing work, your feedback on the methodology/results is most welcome :smiley:

feldmeier commented 1 year ago

HUESTIS is working fine with the latest changes in the QP solver. PRIMAL3 has been working for a while. STADAT1 is still failing, with varying reasons. QGROW22: with the latest changes, it appears to run okay, but makes only very small steps towards optimality. Even after two hours and Millions of iterations it's nowhere near optimal, but still improving without error messages. AUG2DC is working fine, just incredibly slow due to the large initial nullspace dimension. A phase1 algorithm customized for QP might help with this.

Any other observations on Maros Meszaros?

stephane-caron commented 1 year ago

Thanks @feldmeier for taking a look :smiley: These issues appeared when running qpsolvers_benchmark (a benchmark of open source QP solvers to help users compare and select QP solvers). Since you've updated the QP solver we should re-run it for HiGHS :ok_hand: If you are interested feel free to take a look, the benchmark has a single entry script with hopefully enough guidance to be straightforward (all feedback welcome!)

In the current results there were 52 Maros-Meszaros problems not solved by HiGHS, which you can see by running:

$ ./benchmark.py maros_meszaros check_results
[2023-05-22 11:44:27,140] [info] Loading existing results from results/maros_meszaros.csv (results.py:80)
[2023-05-22 11:44:27,148] [info] Check out `results` for the full results data (benchmark.py:233)
[2023-05-22 11:44:27,148] [info] Check out `df` for results as a pandas DataFrame (benchmark.py:235)
Python 3.8.10 (default, Mar 13 2023, 10:26:41) 
Type 'copyright', 'credits' or 'license' for more information
IPython 8.0.1 -- An enhanced Interactive Python. Type '?' for help.

In [1]: df[(df["solver"] == "highs")][(df["found"] == False)][(df["settings"] == "default")]

The output is:

Out[1]: 
       problem solver settings      runtime  found  primal_residual  dual_residual  duality_gap  cost_error
9        AUG2D  highs  default     0.000000  False              inf            inf          inf         NaN
30      AUG2DC  highs  default     0.000000  False              inf            inf          inf         NaN
51    AUG2DCQP  highs  default  1000.145015  False              inf            inf          inf         NaN
72     AUG2DQP  highs  default  1000.120207  False              inf            inf          inf         NaN
156    AUG3DQP  highs  default  1000.011053  False              inf            inf          inf         NaN
177      BOYD1  highs  default  1000.308825  False              inf            inf          inf         NaN
198      BOYD2  highs  default  1000.153471  False              inf            inf          inf         NaN
240   CONT-100  highs  default  1000.058597  False              inf            inf          inf         NaN
261   CONT-101  highs  default  1000.050752  False              inf            inf          inf         NaN
282   CONT-200  highs  default  1163.318025  False              inf            inf          inf         NaN
303   CONT-201  highs  default  1243.499906  False              inf            inf          inf         NaN
324   CONT-300  highs  default     0.000000  False              inf            inf          inf         NaN
555      DTOC3  highs  default  1000.123908  False              inf            inf          inf         NaN
1017  HUES-MOD  highs  default     0.007703  False              inf            inf          inf         NaN
1038   HUESTIS  highs  default     0.008663  False              inf            inf          inf         NaN
1059      KSIP  highs  default     0.003616  False              inf            inf          inf         NaN
1080     LASER  highs  default     0.072293  False              inf            inf          inf         NaN
1374  MOSARQP1  highs  default     1.920444  False              inf            inf          inf         NaN
1605   Q25FV47  highs  default     4.628035  False              inf            inf          inf         NaN
1626  QADLITTL  highs  default  1000.044458  False              inf            inf          inf         NaN
1668    QBANDM  highs  default  1000.004266  False              inf            inf          inf         NaN
1731   QBRANDY  highs  default    13.915137  False              inf            inf          inf         NaN
1773     QE226  highs  default     1.225178  False              inf            inf          inf         NaN
1815  QFFFFF80  highs  default    27.139002  False              inf            inf          inf         NaN
1836  QFORPLAN  highs  default  1000.001666  False              inf            inf          inf         NaN
1857  QGFRDXPN  highs  default     0.221604  False              inf            inf          inf         NaN
1878   QGROW15  highs  default  1000.008930  False              inf            inf          inf         NaN
1899   QGROW22  highs  default     0.104959  False              inf            inf          inf         NaN
1920    QGROW7  highs  default  1000.001176  False              inf            inf          inf         NaN
2046  QPILOTNO  highs  default     4.734245  False              inf            inf          inf         NaN
2130  QSCAGR25  highs  default     0.021207  False              inf            inf          inf         NaN
2172   QSCFXM1  highs  default  1000.009259  False              inf            inf          inf         NaN
2193   QSCFXM2  highs  default     1.910137  False              inf            inf          inf         NaN
2214   QSCFXM3  highs  default  1000.003045  False              inf            inf          inf         NaN
2277    QSCSD1  highs  default  1000.005226  False              inf            inf          inf         NaN
2298    QSCSD6  highs  default     0.266669  False              inf            inf          inf         NaN
2319    QSCSD8  highs  default     0.318647  False              inf            inf          inf         NaN
2340   QSCTAP1  highs  default  1000.001432  False              inf            inf          inf         NaN
2361   QSCTAP2  highs  default  1000.002917  False              inf            inf          inf         NaN
2382   QSCTAP3  highs  default  1000.003934  False              inf            inf          inf         NaN
2424  QSHARE1B  highs  default    69.055565  False              inf            inf          inf         NaN
2445  QSHARE2B  highs  default  1000.014055  False              inf            inf          inf         NaN
2466    QSHELL  highs  default     5.706107  False              inf            inf          inf         NaN
2508  QSHIP04S  highs  default     1.833869  False              inf            inf          inf         NaN
2529  QSHIP08L  highs  default  1000.015269  False              inf            inf          inf         NaN
2571  QSHIP12L  highs  default  1000.022784  False              inf            inf          inf         NaN
2592  QSHIP12S  highs  default     2.168330  False              inf            inf          inf         NaN
2613   QSIERRA  highs  default     0.179179  False              inf            inf          inf         NaN
2634    QSTAIR  highs  default  1000.001925  False              inf            inf          inf         NaN
2697   STADAT1  highs  default     0.000000  False              inf            inf          inf         NaN
2739   STADAT3  highs  default   229.349996  False              inf            inf          inf         NaN
2823      UBH1  highs  default  1000.110255  False              inf            inf          inf         NaN