WIP: scipy callback concept for IPX

[mckib2]

• User callback run every iteration of IPX (excluding crossover)

• TODO

  • pipe more data into scipy_clbk
  • add scipy::clbk_t scipy_clbk as data member if IPM to avoid needless method call-chain passing
  • expose Highs callback data member to consumers

[mckib2]

@mdhaber This is what I meant by modifying our fork of HiGHS. This PR adds a SciPy-specific callback function that will be called every iteration of IPX.

Questions:

• Is this a reasonable approach to deprecating old linprog methods?
• If so, where are reasonable points in the MIP/LP solvers to run callbacks? Every iteration for dual simplex appears to be excessive; after every subproblem for MIP seems reasonable to me

[mdhaber]

This is what I meant by modifying our fork of HiGHS.

Julian said he was willing to take a look. Could you submit a PR to HiGHS?

Is this a reasonable approach to deprecating old linprog methods?

The methods are already deprecated. They can be removed as soon as linprog with method='highs' support callback. I don't think I'm the right person to review the C++, so I can't comment on the approach.

If so, where are reasonable points in the MIP/LP solvers to run callbacks?

They currently run once every iteration of simplex/revised simplex. I think we need to continue to support that. I'd assume we'd do the same for MIP - every iteration of every LP subproblem. If there's something that happens once every MIP subproblem that isn't part of a regular LP iteration, we'd call it once during that, too, presumably with a different phase attribute.

A few other thoughts:

1. Ultimately, I think it would be good if IPX were to use phase=1 for the regular iterations and phase=2 for crossover. But initially, we don't have to support this because there is no crossover in the old linprog solvers. Likewise, we don't (initially) need callback support for MIPs because there is currently not any MIP callback support. If it makes it easier, we do the minimum required to consider the transition backward compatible (aside from the removal of the old methods).

2. If someone uses callback, I don't think they can complain about it slowing things down. They can choose to do nothing on some iterations, even if it's called each iteration.

3. If someone doesn't want to use callback, I can see them complaining that the presence of code written with callback in mind slows things down. I think this is something Julian mentioned being concerned about. If that's the case, I believe it's possible to work around this by branching once at the beginning depending on whether callback is defined.

   • If a callback is not being used, a version of the code runs the way it always used to.
   • If a callback is being used, a new, separate version of the code runs the way it always used to with the exception

   With the appropriate use of function, it seems that this can always be done without much copy-paste, right? For instance, the first version might be a loop that calls a function simplex_iteration each iteration; the second version might be a loop that calls simplex_iteration and callback each iteration. Anyway, the point is to avoid wasting time checking whether callback needs to be called every iteration, if that's considered to be a problem.

4. Will the result object passed into the callback function contain information about the original problem, or a transformed problem? (Preferably the original problem, because it's probably very challenging to document how the transformed problem should be interpreted.)

When you're ready, could you add some minimal Python to linprog so I can take this for a test drive?

[mckib2]

@mdhaber Just created mckib2/scipy#34. I've been testing by running this based on an example in the docs:

from scipy.optimize import linprog
c = [-1, 4]
A = [[-3, 1], [1, 2]]
b = [6, 4]
x0_bounds = (None, None)
x1_bounds = (-3, None)
clbk_fun = lambda res: print("clbk:", res)
linprog(c, A_ub=A, b_ub=b, bounds=[x0_bounds, x1_bounds], callback=clbk_fun, method="highs-ipm", options={"presolve": False})

The only information in the res object currently is nit and it's still only run for IPX (I think I see where to put things in for dual simplex, though). Still figuring out how to get the rest of the information populated in the OptimizeResult object