Mismatch in Row Length and Column Count in Sparse Matrix for TSP Model in Gurobi with R

Question

I'm working on solving a Traveling Salesman Problem (TSP) with time windows using Gurobi in R. My model includes decision variables, subtour elimination constraints, and penalties for early and late arrivals. However, I'm encountering an error when adding constraints to my sparse matrix model$A.

The error message is:

Error in add_row_to_model_A(model, row) : 
  Row length does not match number of columns in model$A

This indicates a mismatch between the row length and the number of columns in model$A. Here is the relevant part of my code:

library(gurobi)
library(Matrix)


data <- data.frame(
  Index = c(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13),
  X = c(0.5, 0.166667, 0.1024, 0.833333, 0.103641, 0.138136, 0.291585, 0.151474, 0.19986, 0.953508, 0.76038, 0.912338, 0.85075, 1),
  Y = c(0.5, 0.833333, 0.166667, 0.333333, 0.839252, 0.833857, 0.146165, 0.160264, 0.181737, 0.397418, 0.371833, 0.287926, 0.212256, 0.339829),
  Predecessor = c(NA, NA, NA, 3, NA, NA, NA, 2, 6, NA, 9, NA, NA, 13),
  Start = c(0, 0, 0, 0.53, 1.3, 1.6, 0, 1.2, 0.33, 1.33, 2.2, 3.43, 3.7, 3.3),
  End = c(6.67, 4.67, 5.33, 2.67, 4.07, 2.93, 4.67, 2.17, 1.5, 6.67, 3.2, 4.73, 5.5, 6.67)
)

# Extract coordinates, predecessors, and time windows
coords <- data.frame(
  index = data$Index,
  x = data$X,
  y = data$Y
)

predecessors <- data$Predecessor
e <- data$Start
l <- data$End

# Create the list of predecessor pairs
predecessor_pairs <- na.omit(data.frame(predecessor = predecessors, index = data$Index))
print("Predecessor pairs:")
print(predecessor_pairs)

# Calculate the distance matrix
distance_matrix <- as.matrix(dist(coords[, c("x", "y")]))
print("Distance matrix:")
print(distance_matrix)

# Define the number of nodes
N <- nrow(coords)

# Define the TSP model using the given formulation
model <- list()

# Add decision variables x_ij
model$varnames <- c()
for (i in 1:N) {
  for (j in 1:N) {
    model$varnames <- c(model$varnames, paste("x", i, j, sep = "_"))
  }
}
model$lb <- rep(0, length(model$varnames))
model$ub <- rep(1, length(model$varnames))
model$vtype <- rep("B", length(model$varnames))

# Add u_i variables for subtour elimination
for (i in 1:N) {
  model$varnames <- c(model$varnames, paste("u", i, sep = "_"))
}
model$lb <- c(model$lb, rep(0, N))
model$ub <- c(model$ub, rep(Inf, N))
model$vtype <- c(model$vtype, rep("C", N))

# Add s_ei and s_li variables for soft time windows
for (i in 1:N) {
  model$varnames <- c(model$varnames, paste("se", i, sep = "_"))
  model$varnames <- c(model$varnames, paste("sl", i, sep = "_"))
}
model$lb <- c(model$lb, rep(0, 2 * N))
model$ub <- c(model$ub, rep(Inf, 2 * N))
model$vtype <- c(model$vtype, rep("C", 2 * N))

# Objective: Minimize the total distance with penalties for early and late arrivals
alpha <- 2
beta <- 4
obj <- c()
for (i in 1:N) {
  for (j in 1:N) {
    obj <- c(obj, distance_matrix[i, j])
  }
}
obj <- c(obj, rep(0, N), rep(alpha, N), rep(beta, N))
model$obj <- obj
model$modelsense <- "min"

# Constraints
model$A <- sparseMatrix(i = NULL, j = NULL, dims = c(0, length(model$varnames)))
model$rhs <- c()
model$sense <- c()

# Helper function to ensure rows are correctly sized
add_row_to_model_A <- function(model, row) {
  print(paste("Adding row of length", length(row), "to matrix with", ncol(model$A), "columns."))
  if (length(row) != ncol(model$A)) {
    stop("Row length does not match number of columns in model$A")
  }
  model$A <- rbind(model$A, row)
}

# Constraint: Sum of x_ij for each i is 1
for (i in 1:N) {
  row <- rep(0, length(model$varnames))
  for (j in 1:N) {
    row[(i-1) * N + j] <- 1
  }
  add_row_to_model_A(model, row)
  model$rhs <- c(model$rhs, 1)
  model$sense <- c(model$sense, "=")
}

# Constraint: Sum of x_ij for each j is 1
for (j in 1:N) {
  row <- rep(0, length(model$varnames))
  for (i in 1:N) {
    row[(i-1) * N + j] <- 1
  }
  add_row_to_model_A(model, row)
  model$rhs <- c(model$rhs, 1)
  model$sense <- c(model$sense, "=")
}

# Constraint: x_ii = 0
for (i in 1:N) {
  row <- rep(0, length(model$varnames))
  row[(i-1) * N + i] <- 1
  add_row_to_model_A(model, row)
  model$rhs <- c(model$rhs, 0)
  model$sense <- c(model$sense, "=")
}

# Calculate M as M ≥ max_{i,j ∈ N} {l_i − e_j + t_ij}
M <- max(sapply(1:N, function(i) sapply(1:N, function(j) l[i] - e[j] + distance_matrix[i, j])))

# Subtour elimination constraints: u_i + t_ij <= u_j + M * (1 - x_ij)
for (i in 1:N) {
  for (j in 1:N) {
    if (i != j) {
      row <- rep(0, length(model$varnames))
      row[length(model$varnames) - 2 * N + i] <- 1  # u_i
      row[length(model$varnames) - 2 * N + j] <- -1  # u_j
      row[(i-1) * N + j] <- M  # x_ij
      add_row_to_model_A(model, row)
      model$rhs <- c(model$rhs, M - distance_matrix[i, j])
      model$sense <- c(model$sense, "<=")
    }
  }
}

# Precedence constraints
for (k in 1:nrow(predecessor_pairs)) {
  pred <- predecessor_pairs$predecessor[k]
  i <- predecessor_pairs$index[k]
  row <- rep(0, length(model$varnames))
  row[length(model$varnames) - 2 * N + pred] <- 1
  row[length(model$varnames) - 2 * N + i] <- -1
  add_row_to_model_A(model, row)
  model$rhs <- c(model$rhs, -1)
  model$sense <- c(model$sense, "<=")
}

# Constraints for soft time windows
for (i in 1:N) {
  # Early arrival penalty
  row <- rep(0, length(model$varnames))
  row[length(model$varnames) - 2 * N + i] <- 1  # u_i
  row[length(model$varnames) - N + i] <- -1  # se_i
  add_row_to_model_A(model, row)
  model$rhs <- c(model$rhs, e[i])
  model$sense <- c(model$sense, ">=")
  
  # Late arrival penalty
  row <- rep(0, length(model$varnames))
  row[length(model$varnames) - 2 * N + i] <- 1  # u_i
  row[length(model$varnames) - N + N + i] <- 1  # sl_i
  add_row_to_model_A(model, row)
  model$rhs <- c(model$rhs, l[i])
  model$sense <- c(model$sense, "<=")
}

# Constraint: u_1 = e_1
row <- rep(0, length(model$varnames))
row[length(model$varnames) - 2 * N + 1] <- 1  # u_1
add_row_to_model_A(model, row)
model$rhs <- c(model$rhs, e[1])
model$sense <- c(model$sense, "=")

# Print model$varnames to ensure consistency
print(paste("Total variables:", length(model$varnames)))

# Solve the model using Gurobi
result <- gurobi(model, list(TimeLimit = 60))

# Extract and print the solution
solution <- result$x
tsp_solution <- matrix(solution[1:(N*N)], nrow = N, ncol = N, byrow = TRUE)
print("TSP Solution Matrix:")
print(tsp_solution)

# Extract the tour
tour <- c()
current_city <- 1
while (length(tour) < N) {
  tour <- c(tour, current_city)
  next_city <- which(tsp_solution[current_city, ] == 1)
  if (length(next_city) == 0 || next_city %in% tour) {
    break
  }
  current_city <- next_city
}
tour <- c(tour, 1)  # Return to the starting city

print("Tour:")
print(tour)

# Calculate the total distance using element-wise multiplication
total_distance <- sum(tsp_solution * distance_matrix)

print(sprintf("Total distance: %.3f", total_distance))

# Calculate the total penalty costs
se_values <- solution[(length(solution) - 2 * N + 1):(length(solution) - N)]
sl_values <- solution[(length(solution) - N + 1):length(solution)]
total_penalty_costs <- sum(alpha * se_values + beta * sl_values)

print(sprintf("Total penalty costs: %.3f", total_penalty_costs))

Question:

Why is there a mismatch between the row length and the number of columns in model$A when adding constraints? How can I ensure that the row lengths always match the number of columns in model$A? Any insights or suggestions to fix this issue would be greatly appreciated!

Debugging with print statements: I added print statements in the add_row_to_model_A function to check the length of the rows being added and the number of columns in model$A. This helped me identify that the mismatch occurs because one of the rows has a length of 238 while the matrix has 237 columns.

Helper function to ensure correct row lengths: I implemented a helper function, add_row_to_model_A, to stop the execution and print a message when there is a length mismatch.

Reviewing the constraint loops: I reviewed the loops where constraints are added to model$A to ensure they generate rows of the correct length. However, I still encounter the mismatch error.