This tutorial was written by Michael Schmitt, a student of the course “Discrete Optimization”.

In this Tutorial, I describe how to solve an ILP using the gomory cut approach implemented in polymake.

This preamble must be integrated in every script file, so that we can use the polytope functionalities of polymake.

use application "polytope";

I introduce a subroutine, that returns the non-integer part of some rational number.

sub f { my $x=shift; my $denom = denominator($x); my $nom = convert_to<Integer>($denom * $x); my $remainder = $nom % $denom; if($remainder < 0){ $remainder += $denom; } return new Rational($remainder, $denom); }

The main Routine of this script gets invoked with the following parameters:

- The Matrix of equalities as an Matrix<Rational>
- The Matrix of inequalites as an Matrix<Rational>
- The Objective function as a list of Rationals.

sub gomory { my $p = new Polytope<Rational>; my $A = shift; my $ineqs = shift; my @obj_list = @_; my $obj = new Vector<Rational>(@obj_list); $p->EQUATIONS = $A; $p->INEQUALITIES = $ineqs;

Since all the inequalities that we add are valid for the integer hull of the polytope, we can conclude that the integer hull is empty, when the examined polytope in any iteration cycle is empty.

if(!$p->FEASIBLE){ print "The integer hull of P is empty.\n"; return; }

$p->LP = new LinearProgram<Rational>(LINEAR_OBJECTIVE=>$obj); my $x = $p->LP->MAXIMAL_VERTEX; print "Optimal solution of the relaxation: ", $x,"\n"; my @B = (); my @N = (); for(my $k = 1; $k < $x->dim(); $k++){ if($x->[$k] > 0){ push @B, $k; }else{ push @N, $k; } }

If we found an integer solution in the linear relaxation, it's the solution of the ILP aswell. If not, we remember the coordinate that is not integer for later steps.

my $i = -1; for(my $k=0; $k<@B; $k++){ if(f($x->[$B[$k]]) != 0){ $i = $k; } } if($i == -1){ print "Found optimal integer solution: ", $x, " -> ", $p->LP->MAXIMAL_VALUE,"\n"; return; }

In this example, I only implemented the computation of a Basis in the non-degenerated case (that is, where every Basis Variables are strictly positive.

my @A_b = (); foreach my $b (@B){ push @A_b, $A->col($b); } my $A_b = transpose(new Matrix<Rational>(@A_b)); my $A_inv = inv($A_b);

Here is where the remembered non-integer variable comes into play, to compute the Gomory cut itself.

my @new_row; my @new_col; for(my $j=0; $j<$A->rows(); $j++){ $new_col[$j] = 0; } $new_row[0] = f($A_inv->row($i)*($A->col(0))); foreach my $j (@N){ $new_row[$j] = (-1)*f((-1)*$A_inv->row($i)*$A->col($j)); } foreach my $j (@B){ $new_row[$j] = 0; } push @new_row, 1; my $new_row_vec = new Vector<Rational>(@new_row); my $new_col_vec = new Vector<Rational>(@new_col); print "New equality: $new_row_vec\n";

The new inequality demands, that the newly introduced slack variable will be positive.

my @new_ineq_col; for(my $j=0; $j<$ineqs->rows(); $j++){ $new_ineq_col[$j] = 0; } my @new_ineq_row = @new_ineq_col; push @new_ineq_row, 0; push @new_ineq_row, 1; my $new_ineq_row_vec = new Vector<Rational>(@new_ineq_row); my $new_ineq_col_vec = new Vector<Rational>(@new_ineq_col);

The two matrices constructed are the equality matrix A and the inequality matrix that restricts the variables to be positive.

In both cases we add a new column on the right side and a new row at the bottom.

my $new_A = new Matrix<Rational>(($A|$new_col_vec) / $new_row_vec); my $new_ineqs = new Matrix<Rational>(($ineqs|$new_ineq_col_vec) / $new_ineq_row_vec); push @obj_list, 0; gomory($new_A, $new_ineqs, @obj_list); }

Finally we create two example matrices and an example objective function, to call the main routine.

my $A = new Matrix<Rational>(<<"."); -14 7 -2 1 0 0 -3 0 1 0 1 0 -3 2 -2 0 0 1 . my $ineqs = new Matrix<Rational>(<<"."); 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 . my @objective = (0,4,-1,0,0,0); gomory($A, $ineqs, @objective);

Except where otherwise noted, content on this wiki is licensed under the following license: CC Attribution-Noncommercial-Share Alike 4.0 International