0

I am trying to find a solution to a traveling salesman problem (TSP) in R. The "TSP" package gives heuristic and often stochastic solutions unless one uses the external Concorde solver. However, on Windows the use of Concorde requires a Cygwin installation. I am working on a server which does not allow me to install Cygwin for security reasons. Is it possible to efficiently find TSP solutions in R without relying on the Concorde solver?

L. Mann
  • 143
  • 6
  • 2
    I solved the TSP problem using this link as guidance. https://dirkschumacher.github.io/ompr/articles/problem-tsp.html – YLC Feb 08 '23 at 22:29
  • @YLC, Did you delete your answer because it requires Cygwin, or ... ? At a glance it looks fine (i.e. not a "link-only" answer) – Ben Bolker Feb 08 '23 at 23:44
  • You could write your own function and use brute force if your data isn't too large. To give you an idea of the sizes-- if you wanted to find the best route through 10 cities, that's 181,440 unique sets (paths you could take). (For n cities, `(n - 1)!/2` routes exist.) – Kat Feb 08 '23 at 23:53
  • @ben-bolker The deleted answer does not require Cygwin. If you believe the deleted answer is useful, I can bring it back. However, I think the link is very informative on its own. I have a leaflet app that uses the code documented in the link. Maybe I should post that as an answer? What do you think? – YLC Feb 09 '23 at 02:38

1 Answers1

0

Use library(ompr) to model the problem as a mixed integer linear program (MILP) using the Miller–Tucker–Zemlin (MTZ) formulation.

Use library(ompr.roi) and library(ROI.plugin.glpk) to solve the problem. GLPK (GNU Linear Programming Kit) solves large-scale MILP problems.

Ref: https://dirkschumacher.github.io/ompr/articles/problem-tsp.html

Here's an application that uses leaflet

library(leaflet)

city.1 <- c(city='Paris',  lat= 48.8566,  lng=2.3522)
city.2 <- c(city='Saumur',  lat= 47.2600, lng=-0.0769)
city.3 <- c(city='Rennes', lat= 48.1147, lng=-1.6794)
city.4 <- c(city='Caen', lat= 49.1800, lng=-0.3700)
city.5 <- c(city='Rouen', lat= 49.4428,  lng=1.0886)
D <- data.frame(rbind(city.1, city.2, city.3, city.4, city.5))
D$lat <- as.numeric(D$lat)
D$lng <- as.numeric(D$lng)

leaflet() %>%
  addTiles() %>% 
  addCircles(data=D, lng = ~lng, lat = ~lat, radius=5000, color='red',
             popup=~city, label=~city)
  
path <- TSP(D)
path


# on the map
  leaflet() %>%
    addTiles() %>% 
    addCircles(data=path, lng = ~lng, lat = ~lat, radius=5000, color='red',
               popup=~city, label=~city) %>%
    addPolylines(data = path, lng = ~lng, lat = ~lat, group = ~group, weight=2)


    
TSP <- function(D) {

  # distance matrix
    library(geosphere)
    D$lat <- as.numeric(D$lat)
    D$lng <- as.numeric(D$lng)
    df <- data.frame(lng= D$lng, lat=D$lat)
    M <- distm(df)


    library(ompr)
    n <- nrow(M) # problem size 
    
    dist_fun <- function(i, j) {
      # a distance extraction function
      vapply(seq_along(i), function(k) M[i[k], j[k]], numeric(1L))
    }
    
    
    model <- MILPModel() %>%
      
      # variable == 1 iff we travel from city i to j
      add_variable(x[i, j], i = 1:n, j = 1:n, type = "integer", lb = 0, ub = 1) %>%
      
      # a helper variable for the MTZ formulation of the tsp
      add_variable(u[i], i = 1:n, lb = 1, ub = n) %>% 
      
      # minimize travel distance
      set_objective(sum_expr(colwise(dist_fun(i, j)) * x[i, j], i = 1:n, j = 1:n), "min") %>%
      
      # you cannot go to the same city
      set_bounds(x[i, i], ub = 0, i = 1:n) %>%
      
      # leave each city
      add_constraint(sum_expr(x[i, j], j = 1:n) == 1, i = 1:n) %>%
      
      # visit each city
      add_constraint(sum_expr(x[i, j], i = 1:n) == 1, j = 1:n) %>%
      
      # ensure no subtours (arc constraints)
      add_constraint(u[i] >= 2, i = 2:n) %>% 
      add_constraint(u[i] - u[j] + 1 <= (n - 1) * (1 - x[i, j]), i = 2:n, j = 2:n)
    
#    model    # builds the model
    
    
  # solve the model
    library(ompr.roi)
    library(ROI.plugin.glpk)
    result <- solve_model(model, with_ROI(solver = "glpk", verbose = TRUE))
    
  # extract the solution
    solution <- get_solution(result, x[i, j])
    path <- solution[which(solution$value > 0),]
    
    odd <- function(x) x%%2 != 0;    even <- function(x) x%%2 == 0
    
    from <- 1
    df <- data.frame()
    for (i in 1:nrow(path)) {
      to <- path$j[path$i == from]
#      print(paste("from", dest$city[from], " to", dest$city[to]))
      temp <- data.frame(city=D$city[from], lat=D$lat[from], lng=D$lng[from], group="A")
      if (even(i)) temp$group[1] <- "B" 
      df <- rbind(df, temp)
      from <- to
    }
    temp <- data.frame(city=D$city[from], lat=D$lat[from], lng=D$lng[from], group="A")
    if (odd(i)) temp$group[1] <- "B" 
    df <- rbind(df, temp)
    
    return(df)

}
  
YLC
  • 104
  • 1
  • 5