# 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))

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)' )