Jupyter Snippet P4M 09

Integer & Linear Programming

An example

Setting up data: cost matrix, demand, supply

F = range(2) # 2 factories
R = range(12) # 12 retailers
C = [[1,2,2,1,3,4,5,7,5,2,3,2],
    [4,5,5,4,1,3,1,2,1,2,4,6]]
demand = [9,4,2,6,4,5,7,8,3,6,9,5]
print('Total demand =',sum(demand))
supply = [ 34, 45]

Total demand = 68

Now we can define a linear program $$\min \sum_{i\in F}\sum_{j\in R} c_{ij} x_{ij}$$ Subject To $$\sum_{r\in R} x_{fr} \le s_r\quad\forall f\in F$$ $$\sum_{f\in F} x_{fr} = d_f\quad\forall r\in R$$ $$x \ge 0$$

from mymip.mycplex import Model

lp = Model()
# define double indexed variables and give them a meaningful names
x = [ [lp.variable('x%dto%d'%(i,j)) for j in R] 
     for i in F ]

lp.min( sum( C[i][j] * x[i][j] for i in F for j in R))
# constraints can be given names too:
lp.SubjectTo({"Supply%d"%f:   sum(x[f][r] for r in R) <= supply[f]   for f in F})
lp.SubjectTo(("Demand%02d"%r, sum(x[f][r] for f in F) == demand[r] ) for r in R)
for f in F: 
    for r in R: x[f][r] >= 0   # all variables non-negative
lp.param["SCRIND"]=1    # set parameter to show CPLEX output 
lp.optimise()

print("The minimum cost is",lp.objective())
for r in R:
    for f in F:
        if x[f][r].x > 0: # amount is not zero
            print("%.1f from F%d to R%02d"%(x[f][r].x,f,r))

Tried aggregator 1 time.
LP Presolve eliminated 0 rows and 4 columns.
Aggregator did 12 substitutions.
Reduced LP has 2 rows, 8 columns, and 16 nonzeros.
Presolve time = 0.01 sec. (0.02 ticks)

Iteration log . . .
Iteration:     1   Dual objective     =           122.000000
The minimum cost is 122.0
9.0 from F0 to R00
4.0 from F0 to R01
2.0 from F0 to R02
6.0 from F0 to R03
4.0 from F1 to R04
5.0 from F1 to R05
7.0 from F1 to R06
8.0 from F1 to R07
3.0 from F1 to R08
6.0 from F1 to R09
8.0 from F0 to R10
1.0 from F1 to R10
5.0 from F0 to R11

To see how the solve sees this problem, try writing it out to file and printing the contents of the file:

lp.write("myfirst.lp")
print(open("myfirst.lp","r").read())

\ENCODING=ISO-8859-1
\Problem name: Model

Minimize
 obj: x0to0 + 2 x0to1 + 2 x0to2 + x0to3 + 3 x0to4 + 4 x0to5 + 5 x0to6 + 7 x0to7
      + 5 x0to8 + 2 x0to9 + 3 x0to10 + 2 x0to11 + 4 x1to0 + 5 x1to1 + 5 x1to2
      + 4 x1to3 + x1to4 + 3 x1to5 + x1to6 + 2 x1to7 + x1to8 + 2 x1to9
      + 4 x1to10 + 6 x1to11
Subject To
 Supply0:  x0to0 + x0to1 + x0to2 + x0to3 + x0to4 + x0to5 + x0to6 + x0to7
           + x0to8 + x0to9 + x0to10 + x0to11 <= 34
 Supply1:  x1to0 + x1to1 + x1to2 + x1to3 + x1to4 + x1to5 + x1to6 + x1to7
           + x1to8 + x1to9 + x1to10 + x1to11 <= 45
 Demand00: x0to0 + x1to0  = 9
 Demand01: x0to1 + x1to1  = 4
 Demand02: x0to2 + x1to2  = 2
 Demand03: x0to3 + x1to3  = 6
 Demand04: x0to4 + x1to4  = 4
 Demand05: x0to5 + x1to5  = 5
 Demand06: x0to6 + x1to6  = 7
 Demand07: x0to7 + x1to7  = 8
 Demand08: x0to8 + x1to8  = 3
 Demand09: x0to9 + x1to9  = 6
 Demand10: x0to10 + x1to10  = 9
 Demand11: x0to11 + x1to11  = 5
End

Advanced Usage

All of the raw CPLEX callable library functions as per the can be accessed if required using C interface (see IBM’s CPLEX documentation) can be accessed by using the cplex object.

For example to identify a minimal conflict set leading to infeasiblity of the problem we can use the CPXrefineconflict and write out the conflict set using CPXclpwrite

from mymip.mycplex import cplex
lp.SubjectTo(x[0][0]+x[1][1] >= 50) # make it infeasible
lp.optimise() # error!
cplex.CPXrefineconflict(lp.Env,lp.LP,0,0) # 2 null pointers
cplex.CPXclpwrite(lp.Env,lp.LP,b"conflict.lp") # note binary (ascii) string
print("#"*10,"Conflict LP","#"*10)
print(open("conflict.lp","r").read())

Iteration log . . .
Iteration:     1   Dual objective     =           204.000000
########## Conflict LP ##########
\ENCODING=ISO-8859-1
\Problem name: Model_conflict

Minimize
 obj:
Subject To
 Demand00: x0to0 + x1to0  = 9
 Demand01: x0to1 + x1to1  = 4
 C15:      x0to0 + x1to1 >= 50
\Sum of equality rows in the conflict:
\ sum_eq: x0to0 + x0to1 + x1to0 + x1to1  = 13
Bounds
      x0to0 Free
      x1to1 Free
End