Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ENH: modularize presolve in linprog #11691 #12510

Closed
wants to merge 38 commits into from
Closed

Conversation

janvle
Copy link
Contributor

@janvle janvle commented Jul 8, 2020

Reference issue

Third (last) step to implement suggestion #11570

What does this implement/fix?

This PR implements _presolve() in a modular way.

Additional information

Andersen & Andersen suggest a large number of presolve possibilities, which were not all implemented in the code (yet). I have not implemented any additional steps, but by breaking up _presolve() into parts, I was able to implement the processing loop A&A suggest.
The attached benchmark results show some large differences (from 10 times slower to 10 times faster). The slower benchmarks have in common that they relate to sparse matrices. The source is the modifying of intermediate sparse matrices (which results in a SparseEfficiencyWarning, see code comments). I don't know if this is sufficiently important to address, and I also don't know how to go about fixing this.
presolve benchmark 01.txt

Because of the new loop, _presolve() is able to reduce problems more thoroughly than before. I have not yet added tests to cover this extra functionality.

'_presolve()' returns an additional presolve_effect variable. Nothing is done with this data yet, but perhaps it is helpful for monitor its functioning.

janvle added 27 commits March 2, 2020 20:57
Added corresponding tests
Revert tox.ini
Needed to pass standard PR checks
Add _presolve subfunctions to detect infeasible constraints, remove variables
which arre fixed by bounds, remove variable which are fixed by row singletons.
Add tests for these functions. Row singleton tests are not complete.
Included asv.conf.json and reverted .gitignore to original.
The code in _presolve() is split into parts, allowing the process to loop until no furter reduction is obtained.
In the redundancy-removal, the initial zero-row test is removed; _presolve() does that now as a final step.
Removed variables, equations, inequalities & redundancy removal result.
Including redundancy-removal result
Resolved conflict in test__linprog_clean_inputs.py and asv.conf.json
Corrected errors in tests mainly. These were written to an earlier version of the presolve routines.
Same reason and same change as the previous commit.
…nfeasibility checks

These calls are still timing bottlenecks.
…esolve_infeasibile_..._constraints()

The benchmarks improve, but having dense intermediate matrices of the same shape as sparse A_eq and/or A_ub defies logic.
Actually a redo. Preferred parameter formats and comments explaining the tests.
Also suppressed efficiency warning in feasibility checks with sparse matrices
Resolved conflicts, small updates to _linprog.py and _linprog_util.py
The _remove_zero()-function is not necessary now that the _presolve() guarantees that there are no zero-rows im the problem matrices.
I modified this function to do a light check, but the warning makes the tests fail.
@janvle
Copy link
Contributor Author

janvle commented Jul 8, 2020

I forgot I modified _remove_zero_rows() in _remove_redundancy.py. The function is not necessary anymore now that _presolve() always removes zero rows from the problem matrices. The modification I made is that the function still tests for zero-rows, and then warns if there are any. The warning causes test failures in test__remove_redundancy.py and test_shgo.py. I removed the warning for now.

@janvle
Copy link
Contributor Author

janvle commented Jul 8, 2020

Getting into details, two test failures left, one that puzzles me:

  • I'm suppressing a SparseEfficiencyWarning; apparently I'm not supposed to do that. I'll change that.
  • Test scipy/optimize/tests/test__shgo.py::test_13_high_sobol fails on a missing file, scipy/optimize/_shgo_lib/sobol_vec.gz. On my local machine this file exists, is found, and the test succeeds. How to repair that?

@janvle
Copy link
Contributor Author

janvle commented Jul 8, 2020

The TravisCI LINT test flags all the cases in test__presolve().py where I tried to neatly align matrix elements -- too bad that has to go.

The presolve step is not consistently slower. In several cases it is also faster. I haven't dug into it deeply, though. Form the names of the tests which perform better, the faster cases may have to do with detecting infeasibility. The slower cases may have to do with looping. Checking for infeasible equality and inequality constraints is relatively heavy, and it happens each loop.

I still have a couple of days before leaving for a 3 week's holiday. What would be most important to check first?

Aligning matrix elements does not fit PEP8 rules.
@mdhaber
Copy link
Contributor

mdhaber commented Jul 8, 2020

neatly align matrix elements

There was talk of relaxing that rule. Let me double check on that.

Update: see gh-12367. I suppose you can wait on this.


# Remove variables which are fixed by singletons in A_eq
if not _presolve_complete(p_stat):
lp, rev, p_stat = _presolve_remove_equation_row_singletons(lp, p_stat, tol)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can we put all of these functions in a list and loop this? Something like:

for fun in fun_list:
    if not _presolve_complete(p_stat):
        lp, rev, p_stat = fun(lp, p_stat, tol)
        if rev is not None:
            revstack.append(rev)

@mdhaber
Copy link
Contributor

mdhaber commented Jul 8, 2020

What would be most important to check first?

First please take a look at image and eliminate the changes you didn't intend to make. Let's get all the CI passing that should be passing.

Next it would be nice to find (in existing tests) or design (probably by working backward from an irreducible problem) a test that is simplified by a few of the functions on each loop and would loop a few times. I'd also like to have a test that is modified by every one of the functions in a single loop.

Third would be to take a look at the speed so we at least know where any increases are coming from.

if complete:
x = np.zeros(lp.c.shape)
else:
A, b, c, _, x0 = _get_Abc(lp, 0)
Copy link
Contributor

@mdhaber mdhaber Jul 8, 2020

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Have you eliminated the need for c0 entirely? If the option 'disp': True will the value of the original objective function at each iteration print correctly? Will the correct value of the objective function be sent to any callback functions? If so, _get_Abc doesn't need to return it, and _linprog_simplex etc... don't need to accept it as an argument; we can remove the idea from the code.

return (_LPProblem(c, A_ub, b_ub, A_eq, b_eq, bounds, x0),
c0, x, revstack, complete, status, message)
presolve_effect[4] = pps[4] - A_eq.shape[0]
return (lp._replace(A_eq=A_eq, b_eq=b_eq),
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A_ub, b_ub don't need to be here?

@janvle
Copy link
Contributor Author

janvle commented Jul 9, 2020

I'll go through most of your comments and questions later.

There are some files I did not intend to change: .gitmodules, doc/scipy-sphinx-theme, doc/source/_static/scipy-mathjax, doc/sphinxext, scipy/_lib/_uarray/scipychanges.patch, scipy/sparse/linalg/dsolve/SuperLU/scipychanges.patch. What do I do to get them out of this pull request?
How do I change the permissions of scipy/optimize/_shgo_lib/sobol_vec.gz? The file is not tracked by git, so I uploaded it manually.

I'm attaching the slow benchmark results:
presolve benchmark 02.txt
I'm not familiar with what regularly comes out of this benchmark, but the list is long and has quite an extreme range. It is worthwhile to explain what is happening; where do I find information about what the individual tests address?

@mdhaber
Copy link
Contributor

mdhaber commented Jul 9, 2020

Re unintended changes:
https://ffeathers.wordpress.com/2020/03/21/removing-an-updated-file-from-a-pr-on-github/
I think you want Case 3.

You can delete the file scipy/optimize/_shgo_lib/sobol_vec.gz if it is not tracked by Git. The error about it was probably caused by changes to the other files.

If you check git status before committing, it gives instructions for adding specific files to a commit and eliminating changes to the rest.

The benchmark problems don't test specific things. They are just a large collection of real-world problems of varying difficulties from Netlib.

@janvle
Copy link
Contributor Author

janvle commented Jul 9, 2020

OK, I'll check that out.

I have been looking at the most extreme cases of the benchmark. The two which report failure, do not fail when I run these benchmarks separately, in Spyder. (I copy the classes, adjust paths to the .npz files , instantiate them and call the timing method.)
optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': False}), 'WOODW') does not fail. The looped presolve implementation is actually about 4 times as fast.
optimize_linprog.Netlib_infeasible.track_netlib_infeasible(('interior-point', {'sparse': False}), 'bgindy') also does not fail, it detects an infeasible problem in about the same time the master branch takes.

The third benchmark, optimize_linprog.Netlib_presolve.time_netlib_presolve('sparse', 'PILOTNOV') takes much more time in the looped presolve since it loops about 8 times, each time eliminating 3 singleton rows. The original code only eliminates 2 singleton rows (I may have missed one?).

Regarding the last couple of cases, e.g. optimize_linprog.Netlib.track_netlib(('revised simplex', {}), 'E226'): from the code I understand that this is just a conditional wrapper for optimize_linprog.Netlib.time_netlib(...), but it seems as if in these last cases this time_netlib() doesn't get called. Could that be?

I found some limited hints about the background of these benchmark data sets: via Google Scholar, David M Gay, "Electronic mail distribution of linear programming test problems". Still it seems a bit strange to me to benchmark code without knowing if these benchmark problems are not biased in some way. But perhaps the large number averages out a lot.

@mdhaber
Copy link
Contributor

mdhaber commented Jul 9, 2020

The two which report failure, do not fail when I run these benchmarks separately, in Spyder.

This is something I've found with ASV, too. Sometimes it reports failure when problems take too long to solve, and maybe that's what's going on here. How long did they take in Spyder?
There is some information on setting timeout limits here. Not sure how to change it for all tests at the same time.

The third benchmark, optimize_linprog.Netlib_presolve.time_netlib_presolve('sparse', 'PILOTNOV') takes much more time in the looped presolve

OK. Yes, I see that most of these large "regressions" are in time_netlib_presolve, and it makes sense that presolve takes longer when it loop. But I need to be convinced that the increase in time is primarily because it is doing a more thorough job - the individual presolve operations should be just as fast as before. You say that

+     3.50±0.01ms       82.4±0.4ms    23.58  optimize_linprog.Netlib_presolve.time_netlib_presolve('sparse', 'PILOTNOV')

is looping 8 times, but it's slowing down presolve by a factor of 24, which suggests that individual iterations are taking much longer.
Then again, the total solution time for PILOTNOV is reduced:

-      2.59±0.01s       2.21±0.01s     0.85  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'PILOTNOV')

so I think the longer presolve is OK. It would just help to know why it is happening.

Regarding the last couple of cases

The track_ benchmarks record the difference in error between the accepted optimal value of the objective and the value returned by linprog. So lines like:

-        1.58e-07         1.08e-09     0.01  optimize_linprog.Netlib.track_netlib(('interior-point', {'sparse': True}), 'RECIPE')

are good because the error was reduced by a factor of 100. Those last few lines:

-            0.96         7.82e-12     0.00  optimize_linprog.Netlib.track_netlib(('revised simplex', {}), 'BRANDY')
-            3.12         1.98e-11     0.00  optimize_linprog.Netlib.track_netlib(('revised simplex', {}), 'E226')
-            0.61         2.13e-16     0.00  optimize_linprog.Netlib.track_netlib(('revised simplex', {}), 'RECIPE')

suggest that revised simplex was basically unable to solve the problem before (large error) and after the improved presolve it solves them with great accuracy.

these benchmark problems are not biased in some way

I think they are intentionally biased toward being numerically challenging, but I'm not sure.


Other thoughts:

There are a lot of big improvements in time_netlib_infeasible. Are these problems being found to be infeasible during presolve?

I'm not concerned by any of the regressions in the track_netlib benchmarks. Some problems get more accurate, some get less accurate, but all the problems that were solved accurately before are still solved accurately, and some problems that were not solved accurately before are now solved accurately. Overall, it looks like solution is more accurate than before.

I'm looking through some of the time_netlib regressions. This one:

+      75.0±0.2ms          389±3ms     5.18  optimize_linprog.Netlib.time_netlib(('revised simplex', {}), 'E226')

is totally fine because it's accompanied by:

-            3.12         1.98e-11     0.00  optimize_linprog.Netlib.track_netlib(('revised simplex', {}), 'E226')

It's OK to take longer now because it actually solves the problem. Same with BRANDY and RECIPE.

So overall these results look OK. I'd just like to better understand what is going on with the large number (45-50) of minor regressions, e.g.:

+         574±7ms          786±4ms     1.37  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'SHIP08L')
+         242±2ms          323±1ms     1.34  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'SHIP12S')
+         199±5ms          261±6ms     1.31  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'GFRD-PNC')
+       144±0.6ms          187±1ms     1.30  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'STANDMPS')
+       150±0.4ms          190±3ms     1.27  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'SHIP04L')
+      95.1±0.4ms          119±1ms     1.25  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'STANDATA')
+         189±2ms          237±3ms     1.25  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'SHIP08S')
+      95.0±0.5ms        118±0.6ms     1.24  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'SHIP04S')
+      4.08±0.02s          4.95±0s     1.21  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'BNL2')
+      1.09±0.01s       1.30±0.01s     1.20  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'SHIP12L')
+      2.21±0.04s       2.62±0.03s     1.19  optimize_linprog.Netlib.time_netlib(('interior-point', {'sparse': True}), 'STOCFOR2')

Is the time difference accounted for by an increase in presolve time?
If so, is that because the presolve goes through multiple iterations, or is it just because presolve is taking longer?
If the presolve is just going through more iterations, but this doesn't make the linear programming algorithm any faster, that's OK. Andersen admits that presolve can make things slower overall in some cases. So we just need to make sure it's as fast as it can be and that we give users the option to control it. What do you think about exposing an option to control the number of presolve iterations?

@@ -783,8 +793,8 @@ def test_unbounded_no_nontrivial_constraints_1(self):
res = linprog(c, A_ub, b_ub, A_eq, b_eq, bounds,
method=self.method, options=self.options)
_assert_unbounded(res)
assert_equal(res.x[-1], np.inf)
assert_equal(res.message[:36], "The problem is (trivially) unbounded")
# assert_equal(res.x[-1], np.inf) # unclear if this should be required
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These were written because the problem should be solved in presolve.

status = 2
message = ("The problem is (trivially) infeasible since one "
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It looks like these messages are gone?

These two directories are empty. Will see if this is all.
Directory is empty. Will see if that is sufficient.
Permissions of the previously uploaded file were wrong (644 rw-r--r-- instead of 755 rwxr-xr-x)
@janvle
Copy link
Contributor Author

janvle commented Jul 10, 2020

Thanks for the explanation, I'm beginning to get a feel for how to look at these benchmarks. I thought all comparisons were times, but some are errors - I didn't now that.

Going through a couple of topics before I'm taking a 3 week's break :-)

I modified my local branch to get rid of the unintentional changes which were present in the PR. I also conformed to PEP8. I succeeded for all but one: not for the sobol_vec.gz file permissions. Uploading changes 755 to 644, and I don't see a way to change that. I see the CircleCI-checks now pass.

The failing benchmarks take 22s (WOODW) and 49s (bgindy), but these numbers may be rather meaningless since the master branch benchmarks differ also: 16s (WOODW) and 47s (bgindy).

An option to limit the number of loops would be OK to me. Perhaps there may be other possibilities, for example restricting presolve to limited functionality. We'll see once we get a clearer view on where extra presolve efforts become a burden.

I'll do these checks when I'm back:

  • if the increased performance in the infeasibility checks is due to presolve
  • what happens to the benchmarks which take slightly longer
  • your review comments

Thanks for reviewing!

@mdhaber
Copy link
Contributor

mdhaber commented Jul 10, 2020

Thanks @janvle. Have a great break.

@tylerjereddy
Copy link
Contributor

@janvle @mdhaber Gentle ping as we approach within ~2 weeks of branching--I can always bump the milestone on this PR if it is not likely to be ready in time--let me know.

@janvle
Copy link
Contributor Author

janvle commented Nov 12, 2020 via email

@tylerjereddy tylerjereddy modified the milestones: 1.6.0, 1.7.0 Nov 13, 2020
@mdhaber mdhaber removed this from the 1.7.0 milestone May 2, 2021
@tupui tupui closed this in #15512 Feb 9, 2022
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement A new feature or improvement scipy.optimize
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

3 participants