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]:

```
# this solution uses NamedArrays, which are arrays indexed over sets for both x and y dimensions.
using JuMP, Clp, NamedArrays
sites = [ 1, 2, 3]
mills = [:A, :B, :C]
cost_per_haul = 4 # don't forget the return trip!
dist = NamedArray( [8 15 50; 10 17 20; 30 26 15], (sites,mills), ("Sites","Mills") )
supply = Dict(zip( sites, [20 30 45] ))
demand = Dict(zip( mills, [30 35 30] ))
m = Model(solver=ClpSolver())
@variable(m, x[sites,mills] >= 0) # x[i,j] is lumber shipped from site i to mill j.
@constraint(m, sup[i in sites], sum(x[i,j] for j in mills) == supply[i] ) # supply constraint
@constraint(m, dem[j in mills], sum(x[i,j] for i in sites) == demand[j] ) # demand constraint
@objective(m, Min, cost_per_haul*sum( x[i,j]*dist[i,j] for i in sites, j in mills ) ) # minimize transportation cost
status = solve(m)
println(status)
# nicely formatted solution
solution = NamedArray( Int[getvalue(x[i,j]) for i in sites, j in mills], (sites,mills), ("Sites","Mills") )
println( solution )
println()
println("Total cost will be \$", getobjectivevalue(m))
```

In [9]:

```
using JuMP, Clp
m = Model(solver=ClpSolver())
# incidence matrix:
A = [ 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 ]
# supply and demand
b = [ 20, 30, 45, -30, -35, -30 ]
# distances
c = [ 8, 15, 50, 10, 17, 20, 30, 26, 15, ]
@variable(m, x[1:9] >= 0)
@constraint(m, A*x .== b)
@objective(m, Min, 4*dot(c,x))
solve(m)
xsol = getvalue(x)
display( reshape(xsol,3,3)' )
```

In [ ]:

```
```