Improve performance of LP sensitivity summary

[odow]

As raised on Discourse: https://discourse.julialang.org/t/performance-of-lp-rhs-perturbation-range/40551

The function for querying LP sensitivity are pretty slow, mainly because they compute the constraint matrix every call:
https://github.com/JuliaOpt/JuMP.jl/blob/7df1796d305aed88ec4f66ea52b8437fa4babe58/src/lp_sensitivity.jl#L190-L192

There is probably plenty of other low-hanging fruit we could improve.

Edit: in fact, it seems the best thing we could do is compute and return a dictionary with the perturbation range for every constraint.

[geoffleyland]

Hi @odow,

sorry to bring this one back to life (I can start a separate issue if that's better), but I've got a case where lp_sensitivity_report fails on the assert in line 78 of lp_sensitivity2.jl - the assertion being @assert size(B, 1) == size(B, 2).

The basis in question has 13 columns and 12 rows. The A matrix has 28 columns and 12 rows, and basic_cols has 13 ones. termination_status(model) is OPTIMAL, and the model was solved with GLPK. The old style "lp_rhs_perturbation_range" works on this model, but the old code doesn't assert that the basis is square. (and it seems reasonable that it should be?)

Any pointers on how I can do to investigate further?

[odow]

Can you send me the MPS file? write_to_file(model, "model.mps")

[odow]

Can't reproduce. It works for me:

model = read_from_file("lp_sensitivity_report_problem.mps")

# lp_sensitivity_report can't handle Interval{} sets. Convert the constraint.
c = first(all_constraints(model, VariableRef, MOI.Interval{Float64}))
c_obj = constraint_object(c)
delete(model, c)
set_lower_bound(c_obj.func, c_obj.set.lower)
set_upper_bound(c_obj.func, c_obj.set.upper)

set_optimizer(model, GLPK.Optimizer)
optimize!(model)
julia> report = lp_sensitivity_report(model)
SensitivityReport(...)

(I haven't attached the MPS file. Not sure if you wanted to share publicly or not.)

If you run these two lines after an optimize! call

JuMP.jl/src/lp_sensitivity2.jl

Lines 75 to 76 in 2baec4c

	std_form = _standard_form_matrix(model)
	basis = _standard_form_basis(model, std_form)

you should be able to check the basis and see if GLPK reports something nosensical.

[geoffleyland]

Thanks @odow!

Yep, I ran those lines just before the call to lp_sensitivity_report, but not just after optimize, and I wonder if this has something to do with us solving the model several times with different objectives. I'll investigate some some more and figure out at which point it goes wonky.