Retrieve the simplex tableau from Gurobi.jl using JuMP

Hello everyone!

I’m trying to retrieve the optimal simplex tableau after solving a linear program using Gurobi. Thanks to Oscar’s help, I learned how to use the GRBBinvRowi function to access individual rows of the inverse basis matrix. Here’s the helper function I’m using:

using JuMP
import Gurobi

"""
    my_GRBBinvRowi(model::Gurobi.Optimizer, i::Int)

Computes a single row of the simplex tableau (row `i` of B⁻¹A).

`i` is the zero-based index of an affine constraint.

See Gurobi documentation:
https://docs.gurobi.com/projects/optimizer/en/current/reference/c/simplex.html#c.GRBBinvRowi
"""
function my_GRBBinvRowi(model::Gurobi.Optimizer, i::Int)
    grb = model
    pNumVars, pNumConstrs = Ref{Cint}(0), Ref{Cint}(0)
    @assert Gurobi.GRBgetintattr(grb, "NumConstrs", pNumConstrs) == 0
    @assert Gurobi.GRBgetintattr(grb, "NumVars", pNumVars) == 0
    len = pNumVars[] + pNumConstrs[]
    ind, val = zeros(Cint, len), zeros(Cdouble, len)
    x = Gurobi.GRBsvec(len, pointer(ind), pointer(val))
    @assert Gurobi.GRBBinvRowi(grb, i, Ref(x)) == 0
    return ind, val
end

I tested the function with this small example:

c = [1, 2]
A = [-1 1; 1 1; 2 -1]
b = [2, 5, 4]

model = direct_model(Gurobi.Optimizer())
@variable(model, x[1:2] >= 0)
@objective(model, Max, c' * x)
@constraint(model, A * x .<= b)
optimize!(model)

grb = backend(model)
ind, val = my_GRBBinvRowi(grb, 0)

The output was:

(Int32[2, 3, 4, 0, 0], [1.5, -0.5, 1.0, -0.5, 1.0])

Now I’m confused by two things:

  1. The index vector contains two zeros (i.e., index 0 appears twice) with different values. What do these duplicate entries mean?
  2. According to Gurobi’s documentation, the result vector may include entries for slack variables, with slack variable for row i having index cols + i. But:
  • How do artificial variables fit in here?
  • What do the two values at index 0 correspond to? Why does this happen? The index 0 should corresponds to the first variable, x[1].

Any insight into interpreting the result of GRBBinvRowi correctly, especially regarding artificial variables and duplicate indices, would be greatly appreciated!

Thanks in advance!

I guess I should have tested before sending to you! Here you go. The vector is sparse, so the missing values are 0.0.

julia> using JuMP

julia> import Gurobi

julia> function my_GRBBinvRowi(model::Gurobi.Optimizer, i::Int)
           pNumVars, pNumConstrs = Ref{Cint}(0), Ref{Cint}(0)
           @assert Gurobi.GRBgetintattr(model, "NumConstrs", pNumConstrs) == 0
           @assert Gurobi.GRBgetintattr(model, "NumVars", pNumVars) == 0
           len = pNumVars[] + pNumConstrs[]
           ind, val = zeros(Cint, len), zeros(Cdouble, len)
           x = Gurobi.GRBsvec(0, pointer(ind), pointer(val))
           pX = Ref(x)
           @assert Gurobi.GRBBinvRowi(model, i, pX) == 0
           modified_x = pX[]
           return ind[1:modified_x.len], val[1:modified_x.len]
       end
my_GRBBinvRowi (generic function with 1 method)

julia> begin
           c = [1, 2]
           A = [-1 1; 1 1; 2 -1]
           b = [2, 5, 4]
           model = direct_model(Gurobi.Optimizer())
           set_silent(model)
           @variable(model, x[1:2] >= 0)
           @objective(model, Max, c' * x)
           @constraint(model, A * x .<= b)
           optimize!(model)
           grb = backend(model)
           ind, val = my_GRBBinvRowi(grb, 0)
       end
Set parameter LicenseID to value 890341
(Int32[2, 3, 4], [1.5, -0.5, 1.0])

I did some similar stuff in SDDP.jl which might help:

1 Like

Thanks, Oscar!

I have a question regarding Gomory mixed-integer cuts, which are typically derived from a linear program in the following standard form:

min  cᵗx  
s.t. Ax ≥ 0  
     x ≥ 0

However, in practice, most real-world problems are not naturally formulated in this form. I noticed that JuMP provides the function:

lp_matrix_data(model::GenericModel{T})

which returns an LPMatrixData{T} struct that captures the LP in the following canonical form:

\begin{aligned}
\min\ & c^\top x + c_0 \\
\text{s.t.}\ & b \le Ax \le d \\
            & x_1 \le x \le x_2
\end{aligned}

This is quite helpful, but it doesn’t match the “standard form” that Gomory cuts typically require, where all constraints are expressed as inequalities with nonnegative slack variables and all variables are constrained to be nonnegative.

My question:

Is there a built-in function or recommended way in JuMP (or MathOptInterface) to automatically convert a general LP model to the standard form needed for Gomory cut generation, i.e., by introducing slack variables and enforcing nonnegativity constraints?

And since we are using direct_model, does Gurobi provide such an operations?

Thanks again!

Is there a built-in function or recommended way in JuMP (or MathOptInterface) to automatically convert a general LP model to the standard form needed for Gomory cut generation, i.e., by introducing slack variables and enforcing nonnegativity constraints?

No

And since we are using direct_model, does Gurobi provide such an operations?

No

It’s a bit complicated, but if you follow this tutorial you can see how to convert into a standard form that’s pretty close: Writing a solver interface · JuMP

It’s do-able, but you’d need to write your Gomory cuts at the MOI level.