1

I have a list of locations and their weights (calculated distances apart) in a matrix. I would like the optimal solution for each location having 3 connections, minimizing total distance.

costs6 <- matrix(c(0,399671,1525211,990914,1689886,1536081,399671,0,1802419,1128519,1964930,1603803,1525211,1802419,0,814942,164677,943489,990914,1128519.4,814942.7,0,953202,565712,1689886,1964930,164677,953202,0, 1004916,1536081,1603803,943489,565712,1004916,0),ncol=6,byrow=TRUE)
plantcap <- rep(3,6)
citydemand <- rep(3,6)
plant.signs <- rep("=",6)
city.signs <- rep("=",6)
lptrans <- lp.transport(costs6,"min",plant.signs,plantcap,city.signs,citydemand)
lptrans$solution
lptrans

This LP solver returns

     [,1] [,2] [,3] [,4] [,5] [,6]
[1,]    3    0    0    0    0    0
[2,]    0    3    0    0    0    0
[3,]    0    0    3    0    0    0
[4,]    0    0    0    3    0    0
[5,]    0    0    0    0    3    0
[6,]    0    0    0    0    0    3

I am wondering if there is a way to max out any Xij at 1, so that the solver will give me three ones in each column/row, rather than one 3 in each column/row? If not, is there another solver I can use to find the solution?

Qwo
  • 15
  • 4

1 Answers1

1

Something like this, setting it up as an LP problem (assuming a symmetric solution matrix)?

library(lpSolve)

costs6 <- matrix(c(0,399671,1525211,990914,1689886,1536081,
                   399671,0,1802419,1128519,1964930,1603803,
                   1525211,1802419,0,814942,164677,943489,
                   990914,1128519.4,814942.7,0,953202,565712,
                   1689886,1964930,164677,953202,0, 1004916,
                   1536081,1603803,943489,565712,1004916,0),ncol=6,byrow=TRUE)
nLoc <- nrow(costs6)
nParams <- sum(1:(nLoc - 1L))
# set up the constraint matrix
# columns are parameters corresponding to the lower triangular of costs6 (by column)
# the first six constraints are for the row/column sums
# the last 15 constraints are for the maximum number of times each path can be used (1)
nConst <- sum(1:nLoc)
mConst <- matrix(0L, nConst, nParams)
mConst[matrix(c(c(combn(1:nLoc, 2)), rep(1:nParams, each = 2)), ncol = 2)] <- 1L
mConst[(nLoc + 1L):nConst,] <- diag(nParams)
lpSol <- lp(
  direction = "min",
  objective.in = unlist(costs6[lower.tri(costs6)]),
  const.mat = mConst,
  const.dir = c(rep("=", nLoc), rep("<=", nParams)),
  const.rhs = c(rep(3L, nLoc), rep(1L, nParams)),
  all.int = TRUE
)
lpSol
#> Success: the objective function is 8688039
# convert the solution to a transport matrix
mSol <- matrix(0, nLoc, nLoc)
mSol[lower.tri(mSol)] <- lpSol$solution
mSol[upper.tri(mSol)] <- t(mSol)[upper.tri(mSol)]
mSol
#>      [,1] [,2] [,3] [,4] [,5] [,6]
#> [1,]    0    1    1    1    0    0
#> [2,]    1    0    0    1    1    0
#> [3,]    1    0    0    0    1    1
#> [4,]    1    1    0    0    0    1
#> [5,]    0    1    1    0    0    1
#> [6,]    0    0    1    1    1    0
jblood94
  • 10,340
  • 1
  • 10
  • 15
  • This is exactly what I am trying to do. Thank you! – Qwo Jul 22 '22 at 12:47
  • This code is working as expected for even-numbered sample sizes, but is freezing up at lpSol when I use an odd-numbered sample (e.g. 13x13 costs matrix). Any idea why this may be? – Qwo Jul 22 '22 at 19:05
  • It's because it's not possible. You're trying to get a symmetric square matrix with zeros along the main diagonal. The sum of all the elements in such a matrix will always be even, but you are trying to force the row/column sums to be 3, yet you have an odd number of rows/columns, and 3 times an odd number is odd. – jblood94 Jul 23 '22 at 02:24
  • As an exercise, start with an empty 5-by-5 matrix. Fill in the main diagonal with zeros, then try to get the row/column sums to all be 3 while keeping the matrix symmetric. – jblood94 Jul 23 '22 at 02:29