Lumber transportation problem (J. Reeb and S. Leavengood)

Millco has three wood mills and is planning three new logging sites. Each mill has a maximum capacity and each logging site can harvest a certain number of truckloads of lumber per day. The cost of a haul is \$2/mile of distance. If distances from logging sites to mills are given below, how should the hauls be routed to minimize hauling costs while meeting all demands?

Logging Site Mill A Mill B Mill C Max loads per day
1 8 15 50 20
2 10 17 20 30
3 30 26 15 45
Mill demand 30 35 30
In [1]:
%load_ext gams_magic
In [2]:
%%capture
%run ~/share/DataTransform.ipynb
from IPython.display import display
In [3]:
%%gams
set sites 'sites', mills 'mills';
parameter dist(sites<,mills<), supply(sites), demand(mills);
scalar cost_per_haul;
In [4]:
# Python Data Input
import pandas as pd
sites = ['1', '2', '3']
mills = ['Mill A','Mill B', 'Mill C']
dist = pd.DataFrame(index=sites, columns=mills, data = [[8, 15, 50], [10, 17, 20], [30, 26, 15]])
supply = list(zip(sites,[20,30,45]))
demand = list(zip(mills,[30,35,30]))
cost_per_haul = 4
display(dist)
dist = gt_from2dim(dist)
Mill A Mill B Mill C
1 8 15 50
2 10 17 20
3 30 26 15
In [5]:
%gams_push dist supply demand cost_per_haul
%gams display dist, supply, demand, cost_per_haul
%gams_lst -e
E x e c u t i o n


----     16 PARAMETER dist  Domain loaded from dist position 2

       Mill A      Mill B      Mill C

1       8.000      15.000      50.000
2      10.000      17.000      20.000
3      30.000      26.000      15.000


----     16 PARAMETER supply  

1 20.000,    2 30.000,    3 45.000


----     16 PARAMETER demand  

Mill A 30.000,    Mill B 35.000,    Mill C 30.000


----     16 PARAMETER cost_per_haul        =        4.000  



In [6]:
%%gams 
# GAMS Model
alias (s,sites), (m,mills);
Positive variable ship(sites,mills);
Variable obj;

Equation defcost;   defcost.. obj =e= cost_per_haul*sum((s,m), ship(s,m)*dist(s,m));
Equation defsupply; defsupply(s).. sum(m, ship(s,m)) =e= supply(s);
Equation defdemand; defdemand(m).. sum(s, ship(s,m)) =e= demand(m);

model millco / all /;
In [7]:
%gams solve millco min obj using lp;
Out[7]:
Solver Status Model Status Objective #equ #var Model Type Solver Solver Time
0 Normal (1) Optimal Global (1) 5760.0 7 10 LP CPLEX 0.002
In [8]:
%gams_pull -d ship obj
display(gt_pivot2d(ship))
print('Total cost will be ' + str(obj['level'][0]))
Mill A Mill B Mill C
1 0.0 20.0 0.0
2 30.0 0.0 0.0
3 0.0 15.0 30.0
Total cost will be 5760.0
In [9]:
%gams_cleanup
import numpy as np
A = np.array([[ 1, 1, 1, 0, 0, 0, 0, 0, 0],
              [ 0, 0, 0, 1, 1, 1, 0, 0, 0],
              [ 0, 0, 0, 0, 0, 0, 1, 1, 1],
              [ 1, 0, 0, 1, 0, 0, 1, 0, 0],
              [ 0, 1, 0, 0, 1, 0, 0, 1, 0],
              [ 0, 0, 1, 0, 0, 1, 0, 0, 1]
             ])
b = np.array([ 20, 30, 45, 30, 35, 30])
c = np.array([ 8, 15, 50, 10, 17, 20, 30, 26, 15])
In [10]:
%%gams
set i,j;
parameter A(i,j), b(i), c(j);

Positive Variable x(j);
Variable obj;

Equation defobj; defobj.. obj =e= 4*sum(j, c(j)*x(j));
Equation e(i);   e(i)..   sum(j, A(i,j)*x(j)) =e= b(i);

model gen / all /;
In [11]:
%gams_push A b c
%gams solve gen min obj using lp;
Out[11]:
Solver Status Model Status Objective #equ #var Model Type Solver Solver Time
0 Normal (1) Optimal Global (1) 5760.0 7 10 LP CPLEX 0.005
In [12]:
# idx interface allows only parameters to be comminicated
%gams Parameter xl(j), objl; xl(j) = x.l(j); objl = obj.l;
%gams_pull -n xl objl
display(xl.reshape(3,3))
print('Total cost will be ' + str(objl))
array([[ 0., 20.,  0.],
       [30.,  0.,  0.],
       [ 0., 15., 30.]])
Total cost will be 5760.0