Lp sensitivity summary #1917
Lp sensitivity summary #1917
Conversation
There was a problem hiding this comment.
Good stuff!
I wonder if converting to standard form can be greatly simplified using a MOIU Model and bridges.
Basically, define new model like
MOIU.@model(MyModel,
(),
(MOI.Interval,),
(),
(),
(MOI.SingleVariable,),
(),
(MOI.Zeros,),
(MOI.VectorAffineFunction)
)and then MOI.copy_to(my_model, JuMP.backend(model)) (you might need to use some automatic bridges). Then this would automatically convert the problem to standard form.
There are quite a few style things: http://www.juliaopt.org/JuMP.jl/latest/style/. In particular, prefer long, explicit names over short abbreviations.
I don't think we're too concerned about performance since this is usually just a one-off.
src/constraints.jl
Outdated
| s.t. Ax == 0 | ||
| L <= x <= U | ||
| """ | ||
| function _std_matrix(model::Model) |
There was a problem hiding this comment.
I agree with Oscar, this should be done by defining the standard form as a MOIU model and letting bridges do the transformation during copy_to
There was a problem hiding this comment.
That would indeed be neat. However, I can't get it to work, I might be missing some thing (not so used with automatic bridges...)
Main issue: Every constraint should get an associated slack variable.
My idea was to copy the model to a model having only interval constraints, (Since you can create an MOI.Interval, from other the other sets.) However, I always get an unsupported constraint error. Moreover, I can't find a bridge that creates intervals (only the ones that splits intervals).
If this was accomplished then the slack bridge would accomplish the rest (I think).
There was a problem hiding this comment.
No there is not any bridge that transforms GreaterThan or LessThan to Interval yet. You can try with
MOIU.@model(MyModel,
(),
(MOI.Interval, MOI.GreaterThan, MOI.LessThan, MOI.EqualTo),
(),
(),
(MOI.SingleVariable,),
(),
(MOI.Zeros,),
(MOI.VectorAffineFunction)
)
Codecov Report
@@ Coverage Diff @@
## master #1917 +/- ##
==========================================
+ Coverage 88.51% 91.28% +2.77%
==========================================
Files 32 32
Lines 4059 4110 +51
==========================================
+ Hits 3593 3752 +159
+ Misses 466 358 -108
Continue to review full report at Codecov.
|
Codecov Report
@@ Coverage Diff @@
## master #1917 +/- ##
==========================================
+ Coverage 88.51% 91.17% +2.65%
==========================================
Files 32 33 +1
Lines 4059 4145 +86
==========================================
+ Hits 3593 3779 +186
+ Misses 466 366 -100
Continue to review full report at Codecov.
|
|
Revised according to most of your comments. Also, (as you noted) I tend to use abbreviations, e.g., var for variable. Thank you for your feedback =) It's very appreciated. Also why do the 32bit version fail to build? |
src/lp_sensitivity.jl
Outdated
| end | ||
| UmX_B = U[basic] - x_B | ||
| LmX_B = L[basic] - x_B | ||
| pos = rho .> 1e-7 |
There was a problem hiding this comment.
I would argue that we should stick to the textbook version and warn if we detect possible numerical instability. Otherwise, these numerical constants need a justification.
There was a problem hiding this comment.
I would claim that this 1e-7 is the textbook version =)
C.f. "Computational Techniques of the Simplex Method" by István Maros section 9.3.4.
Out of context but this regards the values for the feasibility range:
In the latter case a typical values are 1e-8 <= e_f <= 1e-6 if
standard double precision arithmetic with 53 bit mantissa is used.
I could double check what notes he does about other choices or correct the current one (I think he applies the tolerance in the numerator) =)
I,e., I really don't think using 0.0 as tolerance is wise (it might not work even for some simple problems), is one justification a reference to this textbook?
The error is [7] lp_rhs_perturbation_range(::ConstraintRef{Model,MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64},MathOptInterface.EqualTo{Float64}},ScalarShape}) at C:\projects\jump-jl\src\lp_sensitivity.jl:203Because Julia doesn't follow JuMP's style recommendation "A user should see a MethodError only for methods that they called directly.", it's on you to infer which internal assumption of backslash is violated. I would guess that you need to use |
|
Regarding the tolerances used:
Some alternative remedies:
Regarding getting the input (constraint matrix, variable bounds, constraint bounds):
Again, thanks for the feedback =) |
| question how much an objective coefficient or standard form rhs coefficient | ||
| (c.f., [`standard_form_rhs`](@ref)) can change without violating primal or dual | ||
| feasibility of the basic solution. | ||
|
|
There was a problem hiding this comment.
Please add a note that this form of sensitivity analysis appears often in textbooks but has limited use in practice because optimal bases are often degenerate.
There was a problem hiding this comment.
(If you object to this characterization, I'd be interested to hear the argument.)
There was a problem hiding this comment.
True, it depends quite a lot on the instance if these ranges are practically useful, I added a paragraph discussing briefly their limitations.
Adds two functionalities:
perturbation_range_of_feasibility(constraint)perturbation_range_of_optimality(variable)(Targets the issue #1332.)
The choice to present the perturbation (differences) rather than the actual ranges is mainly because the rhs is not necessarily a constant, related to #1890. (It seems sane that these two functions have consistent output).
Are the names too verbose? is e.g.,
range_of_feasibilitybetter?The functions can be optimized and probably simplified, currently the
Lower/UpperBoundRefis not used due to #1892.The perturbation_range_of_feasibility is implemented for affine and variable constraints. However, the interval sets are intentionally not supported (the interpretation of shifting the entire interval seems far-fetched).
Maybe all sensitivity analysis (possibly also
shadow_price) should be moved from constraints.jl, this since they are closely connected (using the same internal functions) and in particularperturbation_range_of_optimalityfeels out of place. To where?Tests to check these functions are also added. The Clp solver (jump-dev/Clp.jl#54) was used to test some larger instances (everything seems fine).
Documentation is intentionally postponed until e.g., the desired output, names, etc. are decided.
A last note, these functions can be heavily optimized if vectorized, i.e., if the ranges of all/many coefficients are desired. So it can be implemented if needed.