As Matt suggested (issue #11515) I've been taking a look at the presolve routine and its reference, the Andersen & Andersen paper.

A&A describes a set of options for trivial simplifications of the problem (some not too trivial, though), numbered (i) - (xvii). The current _presolve function handles just a couple of them:

• zero row in equality constraints (i)
• zero row in inequality constraints (i)
• zero column in both constraints (ii)
• row singleton in equality constraints (v)
• row singleton in inequality constraints (v)
• identical bounds (iv)

A&A also suggest some of these options should be repeated. That currently doesn't happen either. The variables bounds and undo contain the information for reversing the simplification in _postsolve. _postsolve als reverses the problem standardization of _get_Abc using bounds.

I have a couple of considerations and questions with respect to updating _presolve:

1. It would be good to redesign the presolve routine in such a way that we can add options one at a time, in a modular way. Right?
2. A&A warn that some options may be computationally heavy, so it may be good to modify the solver option 'presolve'. Instead of just True or False, we may introduce a numeric level which in- or excludes specific presolve options. Can we do that?
3. Do we know if we can 'match' presolve efforts to an actual solver? What presolve activity is always advantageous? What presolve activities make sense for specific solvers only?
4. I think we should not use bounds as carrier of information for reversing _presolve (in _postsolve). For some A&A options the bounds/undo implementation will not work anyway. But it is also not elegant to use a problem parameter for this.
5. Where in the process should presolve take place? In the current implementation, it is the first step, but the A&A-options relate to the standard form. That means that we have to 'translate' their options, e.g. to include A_ub and b_ub. Would it be better to presolve after the standardization?
6. Some A&A options refer to the dual problem - that is not available at the presolve level. (Is it even used? I did not check the specific linprog solvers). What to do with these A&A-options?
7. The final A&A options (xiv - xvii) address rank deficiency, but in a rather ad-hoc way it seems to me. Linprog currently implements a more thorough procedure to remove redundant equations, but I have also read the comment in _remove_redundancy. What is the proper location to address/check linear dependencies? Immediately before the solver is called? Or can it be done before the standardization in _get_Abc? Is there proof that the standardization never introduces linear dependencies?
8. Analyzing matrix rank (or doing a SVD) is necessary, but it also takes time. Are the solver algorithms able to detect (and perhaps repair) linear dependencies in the course of their process more 'cheaply'? I see quite a bit of literature on this subject, which suggests it is a topic worth exploring, but it is outside of my area of expertise.
9. Am I missing something?

Regarding bounds: I'm confused about the use of None and np.inf at various locations in the code. Do they have different meanings? Perhaps related to the _postsolve process? See 4 above: I would prefer to standardize the use of (-)np.inf for an absent bound.

What do you think of this presolve framework?:

Code each A&A option as a separate function PR, with:

• 1 input: current values of the problem parameters pp = (c, A_eq, b_eq, A_ub, b_ub, u, l, x0)
• 3 outputs: modified values of the problem parameters pp = (c, A_eq, b_eq, A_ub, b_ub, u, l, x0), a function REV to undo this specific modification and a status flag (0 : Modification made, 1 : No modification made, 2 : Problem appears to be infeasible, 3 : Problem appears to be unbounded)

The function REV is pre-parametrized to carry out the restauration. It has:

• 1 input: current values of the criterion vector, solution, and initial guess cp = (c, x, x0)
• 1 output: restored values of the criterion vector, solution, and initial guess cp = (c, x, x0)

The REV-functions can be put on a stack for reversing the modifications, as A&A suggest.

Presolve would now work as follows:

1. input: current values of the problem parameters pp = (c, A_e, b_e, A_f, b_f, u, l, x0)

REVstack = []

3. Perform initial presolve routines

(pp, REV, status) = PR_1(pp)
if status == 0: REVstack.append(REV)
if status == 2 or 3: raise error infeasible or unbound

... and other initial presolve routines in the same way

4. Perform the looped presolve routines

while modifications:
    modifications = False
    (pp, REV, status) = PR_2(pp)
    if status == 0: REVstack.append(REV), modifiations = True
    if status == 2 or 3: raise error infeasible or unbound

... and other presolve routine which must be repeated

5. Perform final presolve routines

(pp, REV, status) = PR_3(pp)
if status == 0: REVstack.append(REV)
if status == 2 or 3: raise error infeasible or unbound

... and other final presolve routines in the same way

6. output: modified values of the problem parameters pp = (c, A_e, b_e, A_f, b_f, u, l, x0)
7. output: REVstack containing the functions to reverse the modifications

Postsolve would work like this:

1. input: criterion vector, solution and initial values of the modified problem cp = (c, x, x0)
2. apply functions in REVstack, last one first

for i in range(len(REVstack),-1,-1):
    cp = REVstack[i](cp)

3. output: criterion vector, solution and initial guess cp = (c, x, x0) in their original context