corrected equation in flp

docs/flp.rst

@@ -66,7 +66,7 @@ This situation and its solution are represented in Figure :ref:`fig-loc`.
 .. _fig-loc:
 
 .. figure:: FIGS/flp.png
-   :scale: 75 %
+   :scale: 50 %
    :align: center
 
    Facility location
@@ -170,7 +170,7 @@ Weak and strong formulations
    Let us consider the facility location problem of the previous section, in a situation where the capacity constraint is not important (any quantity may can be produced at each site).  This is referred to as the *uncapacitated facility location problem*.  One way of modeling this situation is to set the value of :math:`M` in the constraint
 
    .. math::
-      \sum_{i=1}^n x_{ij} \leq M_j y_j & \mbox{ for } j=1,\cdots,m
+      \sum_{i=1}^n x_{ij} \leq M_j y_j  \quad \mbox{ for }  j=1,\cdots,m
 
    as a very large number.  Notice that the formulation is correct even if we omit constraints :math:`x_{ij} \leq d_j y_j, \mbox{ for }  i=1,\cdots,n; j=1,\cdots,m`.  Removing that constraint, the problem may suddenly become very difficult to solve, especially as its size increases; the reason is the *big :math:`M` pitfall*.
 

src/transp_dual.py

@@ -0,0 +1,62 @@
+"""
+transp.py: a model for the transportation problem
+
+Model for solving a transportation problem:
+minimize the total transportation cost for satisfying demand at
+customers, from capacitated facilities.
+
+Data:
+    I - set of customers
+    J - set of facilities
+    c[i,j] - unit transportation cost on arc (i,j)
+    d[i] - demand at node i
+    M[j] - capacity
+
+Copyright (c) by Joao Pedro PEDROSO and Mikio KUBO, 2012
+"""
+
+from pyscipopt import Model, quicksum, multidict, SCIP_PARAMSETTING
+
+I, d = multidict({1:80, 2:270, 3:250, 4:160, 5:180})   # demand
+J, M = multidict({1:500, 2:500, 3:500})                # capacity
+c = {(1,1):4,    (1,2):6,    (1,3):9,  # cost
+     (2,1):5,    (2,2):4,    (2,3):7,
+     (3,1):6,    (3,2):3,    (3,3):4,
+     (4,1):8,    (4,2):5,    (4,3):3,
+     (5,1):10,   (5,2):8,    (5,3):4,
+     }
+
+#Initialize model
+model = Model("transportation")
+model.setPresolve(SCIP_PARAMSETTING.OFF)  # required for retrieving dual information
+
+# Create variables
+x = {}
+
+for i in I:
+    for j in J:
+        x[i,j] = model.addVar(vtype="C", name="x(%s,%s)" % (i,j))
+
+# Demand constraints
+for i in I:
+    model.addCons(quicksum(x[i,j] for j in J if (i,j) in x) == d[i], name="Demand(%s)" % i)
+
+# Capacity constraints
+for j in J:
+    model.addCons(quicksum(x[i,j] for i in I if (i,j) in x) <= M[j], name="Capacity(%s)" % j)
+
+# Objective
+model.setObjective(quicksum(c[i,j]*x[i,j]  for (i,j) in x), "minimize")
+
+model.optimize()
+if model.getStatus() == "optimal":
+    print("Optimal value:", model.getObjVal())
+    EPS = 1.e-6
+    for (i,j) in x:
+        if model.getVal(x[i,j]) > EPS:
+            print("sending quantity %10s from factory %3s to customer %3s" % (model.getVal(x[i,j]),j,i))
+    for c in model.getConss():
+        print("dual of", c.name, ":", model.getDualsolLinear(c))
+else:
+    print("Problem could not be solved to optimality")
+