Algoritmo de Búsqueda Local

Introducción.

Una vez que se tiene una solución para el problema, se puede intentar mejorarla mediante algún procedimiento de búsqueda local. Para cada solución se define un conjunto de soluciones vecinas \(N(s)\). Un procedimiento de Búsqueda Local parte de una solución \(s\), la reemplaza por una solución \(s^∗ \in N(s)\) de menor costo (mejor solución en general) y repite el procedimiento hasta que la solución no pueda ser mejorada. Al terminar, se obtiene una solución localmente óptima respecto a la definición de la vecindad. Para obtener \(s^∗\) puede buscarse la mejor solución de \(N(s)\) (estrategia best improvement) o simplemente tomar la primera solución de \(N(s)\) que mejore el costo (estrategia first improvement).

Usualmente se define \(N(s)\) como las soluciones que pueden obtenerse aplicando a \(s\) alguna regla o procedimiento sencillo llamado movida. Las movidas para el VRP pueden clasificarse en movidas de una ruta y movidas multi-ruta. En las movidas de una ruta los clientes que se visitan en una ruta no varían luego de la aplicación del operador, lo que varía es el orden en que se realizan las visitas. En las movidas multi-ruta, además de cambios en el orden de las visitas suele modificarse el conjunto de clientes visitados en cada ruta.

Operador \(\lambda\) intercambio.

Uno de los operadores más conocidos de búsqueda local para una ruta es el \(\lambda\)-intercambio. Un \(\lambda\)-intercambio consiste en eliminar \(\lambda\) arcos de la solución (\(\lambda >1\)) y reconectar los \(\lambda\) segmentos restantes.

Una solución se dice \(\lambda\)-óptima si no puede ser mejorada utilizando \(\lambda\)-intercambios. Se llama \(\lambda\)-opt a un algoritmo de búsqueda local que utiliza \(\lambda\)-intercambios hasta alcanzar una solución \(\lambda\)-óptima.

En una ruta que visita \(n\) clientes, hay \(\binom{n+1}{\lambda}\) maneras posibles de eliminar \(\lambda\) arcos y, dada una elección de arcos a eliminar, hay \(2^{\lambda-1}(\lambda-1)!\) maneras de reconectar la ruta (incluyendo la reconexión que vuelve a generar la misma ruta). Por lo tanto, la cantidad de \(\lambda\)-intercambios posibles es:

\[2^{\lambda-1}(\lambda-1)! \binom{n+1}{\lambda}\]

Recuerde que:

\[ \binom{n+1}{\lambda} = \frac{(n+1)!}{\lambda!(n+1-\lambda)!} \]

Ejemplo heurística de Inserción + busqueda local con \(\lambda\)-intercambio intraruta para problema VRP en R

Durante las jornadas electorales en el Departamento de Córdoba, la Registraduría Nacional del Estado Civil debe realizar la recolección del material electoral sensible (urnas, formularios, actas y votos físicos) desde todos los municipios hacia el centro de consolidación ubicado en Montería.

Este proceso debe efectuarse inmediatamente después del cierre de las mesas de votación, debido a que:

  • El material es confidencial y de alto valor legal.

  • Debe llegar lo más rápido posible al centro de cómputo departamental.

  • Se busca minimizar riesgos de pérdida, manipulación o retrasos.

  • Los recursos de transporte y seguridad son limitados.

Para esta operación se dispone de una flota de \(12\) vehículos oficiales escoltados, cada uno con capacidad máxima de \(1.5\) toneladas, los cuales parten desde Montería, visitan un subconjunto de municipios para recolectar el material electoral y posteriormente regresan al centro de consolidación.

Cada municipio genera una cantidad estimada de material electoral (expresada en peso o número de cajas), por lo que la asignación de municipios a cada vehículo debe respetar la capacidad máxima de transporte. El peso del material electoral por municipio se encuentra disponible en el siguiente enlace: Peso de material electoral por municipio.

Adicionalmente, se dispone de una matriz de distancias en kilómetros entre municipios: Matriz de distancias en km.

Actualmente, las rutas se planifican de forma empírica, lo que genera:

  • mayores distancias recorridas,

  • incremento en los tiempos de consolidación de resultados,

  • mayores costos de combustible,

  • y mayores riesgos de seguridad por recorridos innecesarios.

En consecuencia, se requiere diseñar un modelo de ruteo de vehículos capacitado (Vehicle Routing Problem, VRP) que permita determinar:

  • qué municipios debe atender cada vehículo,

  • el orden óptimo de visita en cada ruta,

  • y minimizar la distancia total recorrida por la flota.

Datos

library(readxl)
library(dplyr)

#Se lee el archivo de Excel que contiene la matriz de distancias entre municipios
distancias_df<- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/distancias_municipios_cordoba.xlsx")

#Matriz de distancias
municipios <- distancias_df[[1]]  #Extraer la primera columna que contiene los nombres de los municipios
D <- as.matrix(distancias_df[,-1]) #Construir la matriz de distancias eliminando la primera columna
rownames(D) <- municipios #Asignar nombres de filas y columnas a la matriz
colnames(D) <- municipios

#dataframe de demandas
demandas_df <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/demanda.xlsx")

#Matriz de demanda
demandas <- demandas_df$demanda #Crear vector de demandas
names(demandas) <- demandas_df$Municipio #Asignar nombres al vector de demandas para acceder por municipio

#Capacidad del vehículo (toneladas)
Q <- 1.5

#Nodo depósito (punto de origen y regreso de las rutas)
deposito <- "Monteria"

#Tamaño de flota (Número máximo de vehículos disponibles)
K <- 12

Función Inserción Cheapest Insertion Global

cheapest_insertion <- function(D, demandas, Q, K, deposito){
  
  clientes <- setdiff(rownames(D), deposito)  #Se obtiene el conjunto de clientes (todos los nodos excepto el depósito)
  no_asignados <- clientes #Inicialmente todos los clientes están sin asignar
  
  rutas <- list() #Lista que almacenará las rutas
  demanda_ruta <- c() # Vector que almacenará la demanda acumulada por cada ruta
  
  # Ruta inicial
  i <- no_asignados[1] #Se toma el primer cliente no asignado (puede escoger el que quiera o escogerlo aleatorio)
  rutas[[1]] <- c(deposito, i, deposito) #Se crea la primera ruta: depósito -> cliente -> depósito
  demanda_ruta <- c(demandas[i]) #Se inicializa la demanda de esa ruta
  no_asignados <- setdiff(no_asignados, i) #Se elimina ese cliente del conjunto de no asignados
  
  while(length(no_asignados) > 0){ #Mientras existan clientes sin asignar
    
    mejor_delta <- Inf #Inicializar el mejor costo de inserción (delta) como infinito
    mejor_cliente <- NULL #Variables para guardar la mejor decisión encontrada, mejor cliente
    mejor_ruta <- NULL #Variables para guardar la mejor decisión encontrada, mejor ruta
    mejor_pos <- NULL #Variables para guardar la mejor decisión encontrada, mejor posición en la ruta
    
    #Evaluar cada cliente no asignado
    for(k in no_asignados){ 
      
      for(r in 1:length(rutas)){ #Evaluar cada ruta existente
        
        ruta <- rutas[[r]] #Extraer la ruta actual
        
        if(demanda_ruta[r] + demandas[k] <= Q){ #Verificar si el cliente cabe en la ruta (restricción de capacidad)
          
          for(pos in 1:(length(ruta)-1)){  #Evaluar todas las posiciones posibles de inserción en la ruta
            
            i <- ruta[pos] #Nodo anterior en la ruta
            j <- ruta[pos+1] #Nodo siguiente en la ruta
            
            delta <- D[i,k] + D[k,j] - D[i,j] #Calcular el costo incremental de insertar k entre i y j
            
            if(delta < mejor_delta){ #Si esta inserción es mejor (menor costo)
              
              #Actualizar la mejor opción encontrada
              mejor_delta <- delta
              mejor_cliente <- k
              mejor_ruta <- r
              mejor_pos <- pos
              
            }
          }
        }
      }
    }
    
    #Si se encontró una inserción factible
    if(!is.null(mejor_cliente)){
      #Extraer la ruta donde se insertará
      ruta <- rutas[[mejor_ruta]]
      #Actualizar la ruta en la lista
      ruta <- append(ruta, mejor_cliente, after = mejor_pos)
      #Actualizar la ruta en la lista
      rutas[[mejor_ruta]] <- ruta
      #Actualizar la demanda acumulada de la ruta
      demanda_ruta[mejor_ruta] <- demanda_ruta[mejor_ruta] + demandas[mejor_cliente]
      #Eliminar el cliente de los no asignados
      no_asignados <- setdiff(no_asignados, mejor_cliente)
      
    } else { #Si no se pudo insertar en ninguna ruta existente
      
      if(length(rutas) < K){ #Verificar si aún hay vehículos disponibles
        
        nuevo <- no_asignados[1] #Tomar un cliente no asignado
        
        rutas[[length(rutas)+1]] <- c(deposito, nuevo, deposito) #Crear una nueva ruta con ese cliente
        
        demanda_ruta <- c(demanda_ruta, demandas[nuevo]) #Registrar su demanda
        
        no_asignados <- setdiff(no_asignados, nuevo) #Eliminarlo de no asignados
        
      } else { #Si no hay más vehículos, el problema no tiene solución factible
        
        stop("No es posible asignar todos los clientes con la flota disponible")
        
      }
    }
  }
  
  return(rutas)
}

Ejecutar algoritmo

#Se ejecuta la función y se obtienen las rutas resultantes.
rutas <- cheapest_insertion(D, demandas, Q, K, deposito)
rutas

Cantidad de rutas

length(rutas)

Cargas de cada ruta

sapply(rutas, function(r) sum(demandas[r]))

Costo de las rutas

distancia_ruta <- function(ruta, D){
  
  distancia <- 0
  
  for(i in 1:(length(ruta)-1)){
    
    distancia <- distancia + D[ruta[i], ruta[i+1]]
    
  }
  
  return(distancia)
}

distancias <- sapply(rutas, distancia_ruta, D)
distancia_total <- sum(distancias)
distancia_total

Resumen de resultados

data.frame(
  Ruta = 1:length(rutas),
  Distancia_km = sapply(rutas, distancia_ruta, D),
  Carga = sapply(rutas, function(r) sum(demandas[r]))
)

Mapa

#Datos
library(readxl)
library(dplyr)
library(leaflet)

#Se lee el archivo de coordenadas geográficas de los municipios.
coords_df <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/coordenadas.xlsx")
coords_df$Municipio <- toupper(trimws(coords_df$Municipio)) #Estandarizan los nombres de los municipios
municipios <- toupper(trimws(municipios)) #Estandarizan los nombres del vector municipios


#Preparación de rutas 

#lapply aplica una función a cada ruta para convertirla en coordenadas.
rutas_coords <- lapply(rutas, function(ruta){ #contiene las rutas generadas por Clarke & Wright.
  
  ruta_municipios <- toupper(ruta) #Convierte los nombres de los municipios de la ruta a mayúsculas para garantizar coincidencia con coords_df.
  
  coords_df %>%
    slice(match(ruta_municipios, Municipio)) #Encuentra la posición de cada municipio de la ruta dentro del dataframe de coordenadas. Slece Selecciona esas filas en ese mismo orden.
  
})

#Paleta de colores
colores <- colorFactor(rainbow(length(rutas)), domain = 1:length(rutas)) #Se genera una función de colores para diferenciar rutas. rainbow() crea colores distintos, length(rutas) define cuántos colores se necesitan. 

#Crear mapa: Se crea un mapa interactivo vacío.
mapa <- leaflet() %>%
  addTiles()

#Dibujar cada ruta
for(i in 1:length(rutas_coords)){ #Se recorre cada ruta generada por el algoritmo.
  
  ruta_df <- rutas_coords[[i]] #Se obtiene el dataframe de coordenadas de la ruta i.
  
  mapa <- mapa %>%
    addPolylines(
      lng = ruta_df$longitud,
      lat = ruta_df$latitud,
      color = colores(i),
      weight = 4
    ) %>%
    
    addCircleMarkers(
      lng = ruta_df$longitud,
      lat = ruta_df$latitud,
      radius = 5,
      color = colores(i),
      popup = paste("Vehículo", i, "<br>", ruta_df$Municipio)
    )
}

mapa

Búsqueda local con \(\lambda=2\)

Una vez construida una solución inicial mediante la heurística de inserción, se procede a mejorar cada ruta utilizando el operador \(\lambda\)-intercambio con \(\lambda = 2\), conocido como 2-opt.

Este operador consiste en eliminar dos arcos de una ruta y reconectar los segmentos resultantes invirtiendo el orden de visita de un subconjunto de nodos. La implementación empleada sigue una estrategia de tipo first improvement, en la cual se acepta el primer movimiento que reduzca el costo de la ruta.

Es importante destacar que este procedimiento actúa únicamente sobre cada ruta de manera independiente (intra-ruta), por lo que no modifica la asignación de clientes entre vehículos.

# Función 2-opt
two_opt <- function(ruta, D, max_iter = 10000000){
  
  mejoro <- TRUE
  iter <- 0
  
  while(mejoro && iter < max_iter){
    
    iter <- iter + 1
    mejoro <- FALSE
    n <- length(ruta)
    
    # Si la ruta es muy corta, no hay nada que optimizar
    if(n <= 4) break
    
    for(i in 2:(n-3)){
      
      for(j in seq(i+2, n-2)){  #evita j = i+1 automáticamente
        
        if(j+1 > n) next
        
        a <- ruta[i]
        b <- ruta[i+1]
        c <- ruta[j]
        d <- ruta[j+1]
        
        # Validación contra la matriz
        if(!(a %in% rownames(D) && 
             b %in% rownames(D) && 
             c %in% rownames(D) && 
             d %in% rownames(D))) next
        
        delta <- D[a,c] + D[b,d] - D[a,b] - D[c,d]
        
        if(delta < -1e-6){
          
          ruta[(i+1):j] <- rev(ruta[(i+1):j])
          mejoro <- TRUE
          break
        }
      }
      
      if(mejoro) break
    }
  }
  
  return(ruta)
}

Función aplicación de función 2-opt a todas las rutas

busqueda_local_2opt <- function(rutas, D){
  
  rutas_mejoradas <- lapply(rutas, function(ruta){
    
    # Se mejora cada ruta de manera independiente
    two_opt(ruta, D)
    
  })
  
  return(rutas_mejoradas)
}

Ejecución 2-opt

# Aplicar búsqueda local con λ = 2
rutas_2opt <- busqueda_local_2opt(rutas, D)

Evaluación de la mejora

# Comparación de distancias antes y después

distancias_original <- sapply(rutas, distancia_ruta, D)
distancias_mejoradas <- sapply(rutas_2opt, distancia_ruta, D)

data.frame(
  Ruta = 1:length(rutas),
  Distancia_Original = distancias_original,
  Distancia_Mejorada = distancias_mejoradas,
  Mejora = distancias_original - distancias_mejoradas
)

# Distancia total
sum(distancias_original)
sum(distancias_mejoradas)

Búsqueda local con \(\lambda=3\)

# Función 3-opt (versión práctica)
three_opt <- function(ruta, D, max_iter = 50){
  
  mejoro <- TRUE
  iter <- 0
  
  while(mejoro && iter < max_iter){
    
    iter <- iter + 1
    mejoro <- FALSE
    n <- length(ruta)
    
    # Rutas muy cortas no se pueden mejorar
    if(n <= 5) break
    
    for(i in 2:(n-4)){
      
      for(j in (i+1):(n-3)){
        
        for(k in (j+1):(n-2)){
          
          a <- ruta[i]
          b <- ruta[i+1]
          c <- ruta[j]
          d <- ruta[j+1]
          e <- ruta[k]
          f <- ruta[k+1]
          
          # Validación
          if(!(a %in% rownames(D) && b %in% rownames(D) &&
               c %in% rownames(D) && d %in% rownames(D) &&
               e %in% rownames(D) && f %in% rownames(D))) next
          
          # Costo actual
          costo_actual <- D[a,b] + D[c,d] + D[e,f]
          
          # ==== OPCIONES DE RECONEXIÓN ====
          
          # Opción 1: invertir segmento (b..c)
          nueva1 <- c(ruta[1:i], rev(ruta[(i+1):j]), ruta[(j+1):n])
          costo1 <- distancia_ruta(nueva1, D)
          
          # Opción 2: invertir segmento (d..e)
          nueva2 <- c(ruta[1:j], rev(ruta[(j+1):k]), ruta[(k+1):n])
          costo2 <- distancia_ruta(nueva2, D)
          
          # Opción 3: doble inversión
          nueva3 <- c(ruta[1:i], rev(ruta[(i+1):j]), 
                      rev(ruta[(j+1):k]), ruta[(k+1):n])
          costo3 <- distancia_ruta(nueva3, D)
          
          # Seleccionar mejor movimiento
          costos <- c(costo1, costo2, costo3)
          mejor_costo <- min(costos)
          
          if(mejor_costo < distancia_ruta(ruta, D) - 1e-6){
            
            if(mejor_costo == costo1) ruta <- nueva1
            if(mejor_costo == costo2) ruta <- nueva2
            if(mejor_costo == costo3) ruta <- nueva3
            
            mejoro <- TRUE
            break
          }
        }
        if(mejoro) break
      }
      if(mejoro) break
    }
  }
  
  return(ruta)
}

Aplicación a todas las rutas

busqueda_local_3opt <- function(rutas, D){
  
  rutas_mejoradas <- lapply(rutas, function(ruta){
    three_opt(ruta, D)
  })
  
  return(rutas_mejoradas)
}

Ejecución

rutas_3opt <- busqueda_local_3opt(rutas, D)

Evaluación de la mejora

# Comparación de distancias antes y después

distancias_original <- sapply(rutas, distancia_ruta, D)
distancias_mejoradas <- sapply(rutas_3opt, distancia_ruta, D)

data.frame(
  Ruta = 1:length(rutas),
  Distancia_Original = distancias_original,
  Distancia_Mejorada = distancias_mejoradas,
  Mejora = distancias_original - distancias_mejoradas
)

# Distancia total
sum(distancias_original)
sum(distancias_mejoradas)

SWAP inter-rutas

La idea de la conformación de la vecindad en este tipo de búsqueda, corresponde a intercambiar nodos entre rutas de la solución \(s\), siempre y cuando no se violen las restricciones, en este caso, la capacidad de los vehículos pertenecientes a la flota.

Ejemplo heurística de Inserción + busqueda local con SWAP inter-ruta para problema VRP en R

Durante las jornadas electorales en el Departamento de Córdoba, la Registraduría Nacional del Estado Civil debe realizar la recolección del material electoral sensible (urnas, formularios, actas y votos físicos) desde todos los municipios hacia el centro de consolidación ubicado en Montería.

Este proceso debe efectuarse inmediatamente después del cierre de las mesas de votación, debido a que:

  • El material es confidencial y de alto valor legal.

  • Debe llegar lo más rápido posible al centro de cómputo departamental.

  • Se busca minimizar riesgos de pérdida, manipulación o retrasos.

  • Los recursos de transporte y seguridad son limitados.

Para esta operación se dispone de una flota de \(12\) vehículos oficiales escoltados, cada uno con capacidad máxima de \(1.5\) toneladas, los cuales parten desde Montería, visitan un subconjunto de municipios para recolectar el material electoral y posteriormente regresan al centro de consolidación.

Cada municipio genera una cantidad estimada de material electoral (expresada en peso o número de cajas), por lo que la asignación de municipios a cada vehículo debe respetar la capacidad máxima de transporte. El peso del material electoral por municipio se encuentra disponible en el siguiente enlace: Peso de material electoral por municipio.

Adicionalmente, se dispone de una matriz de distancias en kilómetros entre municipios: Matriz de distancias en km.

Actualmente, las rutas se planifican de forma empírica, lo que genera:

  • mayores distancias recorridas,

  • incremento en los tiempos de consolidación de resultados,

  • mayores costos de combustible,

  • y mayores riesgos de seguridad por recorridos innecesarios.

En consecuencia, se requiere diseñar un modelo de ruteo de vehículos capacitado (Vehicle Routing Problem, VRP) que permita determinar:

  • qué municipios debe atender cada vehículo,

  • el orden óptimo de visita en cada ruta,

  • y minimizar la distancia total recorrida por la flota.

Datos

library(readxl)
library(dplyr)

#Se lee el archivo de Excel que contiene la matriz de distancias entre municipios
distancias_df<- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/distancias_municipios_cordoba.xlsx")

#Matriz de distancias
municipios <- distancias_df[[1]]  #Extraer la primera columna que contiene los nombres de los municipios
D <- as.matrix(distancias_df[,-1]) #Construir la matriz de distancias eliminando la primera columna
rownames(D) <- municipios #Asignar nombres de filas y columnas a la matriz
colnames(D) <- municipios

#dataframe de demandas
demandas_df <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/demanda.xlsx")

#Matriz de demanda
demandas <- demandas_df$demanda #Crear vector de demandas
names(demandas) <- demandas_df$Municipio #Asignar nombres al vector de demandas para acceder por municipio

#Capacidad del vehículo (toneladas)
Q <- 1.5

#Nodo depósito (punto de origen y regreso de las rutas)
deposito <- "Monteria"

#Tamaño de flota (Número máximo de vehículos disponibles)
K <- 12

Función Inserción Cheapest Insertion Global

cheapest_insertion <- function(D, demandas, Q, K, deposito){
  
  clientes <- setdiff(rownames(D), deposito)  #Se obtiene el conjunto de clientes (todos los nodos excepto el depósito)
  no_asignados <- clientes #Inicialmente todos los clientes están sin asignar
  
  rutas <- list() #Lista que almacenará las rutas
  demanda_ruta <- c() # Vector que almacenará la demanda acumulada por cada ruta
  
  # Ruta inicial
  i <- no_asignados[1] #Se toma el primer cliente no asignado (puede escoger el que quiera o escogerlo aleatorio)
  rutas[[1]] <- c(deposito, i, deposito) #Se crea la primera ruta: depósito -> cliente -> depósito
  demanda_ruta <- c(demandas[i]) #Se inicializa la demanda de esa ruta
  no_asignados <- setdiff(no_asignados, i) #Se elimina ese cliente del conjunto de no asignados
  
  while(length(no_asignados) > 0){ #Mientras existan clientes sin asignar
    
    mejor_delta <- Inf #Inicializar el mejor costo de inserción (delta) como infinito
    mejor_cliente <- NULL #Variables para guardar la mejor decisión encontrada, mejor cliente
    mejor_ruta <- NULL #Variables para guardar la mejor decisión encontrada, mejor ruta
    mejor_pos <- NULL #Variables para guardar la mejor decisión encontrada, mejor posición en la ruta
    
    #Evaluar cada cliente no asignado
    for(k in no_asignados){ 
      
      for(r in 1:length(rutas)){ #Evaluar cada ruta existente
        
        ruta <- rutas[[r]] #Extraer la ruta actual
        
        if(demanda_ruta[r] + demandas[k] <= Q){ #Verificar si el cliente cabe en la ruta (restricción de capacidad)
          
          for(pos in 1:(length(ruta)-1)){  #Evaluar todas las posiciones posibles de inserción en la ruta
            
            i <- ruta[pos] #Nodo anterior en la ruta
            j <- ruta[pos+1] #Nodo siguiente en la ruta
            
            delta <- D[i,k] + D[k,j] - D[i,j] #Calcular el costo incremental de insertar k entre i y j
            
            if(delta < mejor_delta){ #Si esta inserción es mejor (menor costo)
              
              #Actualizar la mejor opción encontrada
              mejor_delta <- delta
              mejor_cliente <- k
              mejor_ruta <- r
              mejor_pos <- pos
              
            }
          }
        }
      }
    }
    
    #Si se encontró una inserción factible
    if(!is.null(mejor_cliente)){
      #Extraer la ruta donde se insertará
      ruta <- rutas[[mejor_ruta]]
      #Actualizar la ruta en la lista
      ruta <- append(ruta, mejor_cliente, after = mejor_pos)
      #Actualizar la ruta en la lista
      rutas[[mejor_ruta]] <- ruta
      #Actualizar la demanda acumulada de la ruta
      demanda_ruta[mejor_ruta] <- demanda_ruta[mejor_ruta] + demandas[mejor_cliente]
      #Eliminar el cliente de los no asignados
      no_asignados <- setdiff(no_asignados, mejor_cliente)
      
    } else { #Si no se pudo insertar en ninguna ruta existente
      
      if(length(rutas) < K){ #Verificar si aún hay vehículos disponibles
        
        nuevo <- no_asignados[1] #Tomar un cliente no asignado
        
        rutas[[length(rutas)+1]] <- c(deposito, nuevo, deposito) #Crear una nueva ruta con ese cliente
        
        demanda_ruta <- c(demanda_ruta, demandas[nuevo]) #Registrar su demanda
        
        no_asignados <- setdiff(no_asignados, nuevo) #Eliminarlo de no asignados
        
      } else { #Si no hay más vehículos, el problema no tiene solución factible
        
        stop("No es posible asignar todos los clientes con la flota disponible")
        
      }
    }
  }
  
  return(rutas)
}

Ejecutar algoritmo

#Se ejecuta la función y se obtienen las rutas resultantes.
rutas <- cheapest_insertion(D, demandas, Q, K, deposito)
rutas

Cantidad de rutas

length(rutas)

Cargas de cada ruta

sapply(rutas, function(r) sum(demandas[r]))

Costo de las rutas

distancia_ruta <- function(ruta, D){
  
  distancia <- 0
  
  for(i in 1:(length(ruta)-1)){
    
    distancia <- distancia + D[ruta[i], ruta[i+1]]
    
  }
  
  return(distancia)
}

distancias <- sapply(rutas, distancia_ruta, D)
distancia_total <- sum(distancias)
distancia_total

Resumen de resultados

data.frame(
  Ruta = 1:length(rutas),
  Distancia_km = sapply(rutas, distancia_ruta, D),
  Carga = sapply(rutas, function(r) sum(demandas[r]))
)

Mapa

#Datos
library(readxl)
library(dplyr)
library(leaflet)

#Se lee el archivo de coordenadas geográficas de los municipios.
coords_df <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/coordenadas.xlsx")
coords_df$Municipio <- toupper(trimws(coords_df$Municipio)) #Estandarizan los nombres de los municipios
municipios <- toupper(trimws(municipios)) #Estandarizan los nombres del vector municipios


#Preparación de rutas 

#lapply aplica una función a cada ruta para convertirla en coordenadas.
rutas_coords <- lapply(rutas, function(ruta){ #contiene las rutas generadas por Clarke & Wright.
  
  ruta_municipios <- toupper(ruta) #Convierte los nombres de los municipios de la ruta a mayúsculas para garantizar coincidencia con coords_df.
  
  coords_df %>%
    slice(match(ruta_municipios, Municipio)) #Encuentra la posición de cada municipio de la ruta dentro del dataframe de coordenadas. Slece Selecciona esas filas en ese mismo orden.
  
})

#Paleta de colores
colores <- colorFactor(rainbow(length(rutas)), domain = 1:length(rutas)) #Se genera una función de colores para diferenciar rutas. rainbow() crea colores distintos, length(rutas) define cuántos colores se necesitan. 

#Crear mapa: Se crea un mapa interactivo vacío.
mapa <- leaflet() %>%
  addTiles()

#Dibujar cada ruta
for(i in 1:length(rutas_coords)){ #Se recorre cada ruta generada por el algoritmo.
  
  ruta_df <- rutas_coords[[i]] #Se obtiene el dataframe de coordenadas de la ruta i.
  
  mapa <- mapa %>%
    addPolylines(
      lng = ruta_df$longitud,
      lat = ruta_df$latitud,
      color = colores(i),
      weight = 4
    ) %>%
    
    addCircleMarkers(
      lng = ruta_df$longitud,
      lat = ruta_df$latitud,
      radius = 5,
      color = colores(i),
      popup = paste("Vehículo", i, "<br>", ruta_df$Municipio)
    )
}

mapa

Función SWAP inter - ruta

#Función swap inter-ruta (first improvement)
swap_interruta <- function(rutas, D, demandas, Q, max_iter = 100){
  
  rutas <- lapply(rutas, function(r) r)  # copia
  mejoro <- TRUE
  iter <- 0
  
  # Función auxiliar: carga de una ruta
  carga_ruta <- function(r){
    sum(demandas[r])
  }
  
  while(mejoro && iter < max_iter){
    
    iter <- iter + 1
    mejoro <- FALSE
    
    for(r1 in 1:(length(rutas)-1)){
      
      for(r2 in (r1+1):length(rutas)){
        
        ruta1 <- rutas[[r1]]
        ruta2 <- rutas[[r2]]
        
        # excluir depósito (primer y último nodo)
        for(i in 2:(length(ruta1)-1)){
          
          for(j in 2:(length(ruta2)-1)){
            
            cliente1 <- ruta1[i]
            cliente2 <- ruta2[j]
            
            # nuevas rutas (swap)
            nueva1 <- ruta1
            nueva2 <- ruta2
            
            nueva1[i] <- cliente2
            nueva2[j] <- cliente1
            
            # verificar capacidad
            if(carga_ruta(nueva1) > Q || carga_ruta(nueva2) > Q) next
            
            # calcular costos
            costo_actual <- distancia_ruta(ruta1, D) + distancia_ruta(ruta2, D)
            costo_nuevo  <- distancia_ruta(nueva1, D) + distancia_ruta(nueva2, D)
            
            # criterio de mejora
            if(costo_nuevo < costo_actual - 1e-6){
              
              rutas[[r1]] <- nueva1
              rutas[[r2]] <- nueva2
              
              mejoro <- TRUE
              break
            }
          }
          
          if(mejoro) break
        }
        
        if(mejoro) break
      }
      
      if(mejoro) break
    }
  }
  
  return(rutas)
}

Ejecución

rutas_swap <- swap_interruta(rutas, D, demandas, Q)

Evaluación

distancias_swap <- sapply(rutas_swap, distancia_ruta, D)

data.frame(
  Ruta = 1:length(rutas),
  Original = sapply(rutas, distancia_ruta, D),
  Swap = distancias_swap,
  Mejora = sapply(rutas, distancia_ruta, D) - distancias_swap
)

sum(sapply(rutas, distancia_ruta, D))
sum(distancias_swap)

Cantidad de rutas

length(rutas_swap)

Cargas de cada ruta

sapply(rutas_swap, function(r) sum(demandas[r]))

Relocation

La idea es:

  • Un cliente \(i\) de una ruta \(R_1\)
  • Se remueve \(i\) de \(R_1\)
  • Se inserta \(i\) en la ruta \(R_2\) en la mejor posición posible.

Repetando capacidad de los vehículos.

Ejemplo heurística de Inserción + busqueda local con Relocation inter-ruta para problema VRP en R

Durante las jornadas electorales en el Departamento de Córdoba, la Registraduría Nacional del Estado Civil debe realizar la recolección del material electoral sensible (urnas, formularios, actas y votos físicos) desde todos los municipios hacia el centro de consolidación ubicado en Montería.

Este proceso debe efectuarse inmediatamente después del cierre de las mesas de votación, debido a que:

  • El material es confidencial y de alto valor legal.

  • Debe llegar lo más rápido posible al centro de cómputo departamental.

  • Se busca minimizar riesgos de pérdida, manipulación o retrasos.

  • Los recursos de transporte y seguridad son limitados.

Para esta operación se dispone de una flota de \(12\) vehículos oficiales escoltados, cada uno con capacidad máxima de \(1.5\) toneladas, los cuales parten desde Montería, visitan un subconjunto de municipios para recolectar el material electoral y posteriormente regresan al centro de consolidación.

Cada municipio genera una cantidad estimada de material electoral (expresada en peso o número de cajas), por lo que la asignación de municipios a cada vehículo debe respetar la capacidad máxima de transporte. El peso del material electoral por municipio se encuentra disponible en el siguiente enlace: Peso de material electoral por municipio.

Adicionalmente, se dispone de una matriz de distancias en kilómetros entre municipios: Matriz de distancias en km.

Actualmente, las rutas se planifican de forma empírica, lo que genera:

  • mayores distancias recorridas,

  • incremento en los tiempos de consolidación de resultados,

  • mayores costos de combustible,

  • y mayores riesgos de seguridad por recorridos innecesarios.

En consecuencia, se requiere diseñar un modelo de ruteo de vehículos capacitado (Vehicle Routing Problem, VRP) que permita determinar:

  • qué municipios debe atender cada vehículo,

  • el orden óptimo de visita en cada ruta,

  • y minimizar la distancia total recorrida por la flota.

Datos

library(readxl)
library(dplyr)

#Se lee el archivo de Excel que contiene la matriz de distancias entre municipios
distancias_df<- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/distancias_municipios_cordoba.xlsx")

#Matriz de distancias
municipios <- distancias_df[[1]]  #Extraer la primera columna que contiene los nombres de los municipios
D <- as.matrix(distancias_df[,-1]) #Construir la matriz de distancias eliminando la primera columna
rownames(D) <- municipios #Asignar nombres de filas y columnas a la matriz
colnames(D) <- municipios

#dataframe de demandas
demandas_df <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/demanda.xlsx")

#Matriz de demanda
demandas <- demandas_df$demanda #Crear vector de demandas
names(demandas) <- demandas_df$Municipio #Asignar nombres al vector de demandas para acceder por municipio

#Capacidad del vehículo (toneladas)
Q <- 1.5

#Nodo depósito (punto de origen y regreso de las rutas)
deposito <- "Monteria"

#Tamaño de flota (Número máximo de vehículos disponibles)
K <- 12

Función Inserción Cheapest Insertion Global

cheapest_insertion <- function(D, demandas, Q, K, deposito){
  
  clientes <- setdiff(rownames(D), deposito)  #Se obtiene el conjunto de clientes (todos los nodos excepto el depósito)
  no_asignados <- clientes #Inicialmente todos los clientes están sin asignar
  
  rutas <- list() #Lista que almacenará las rutas
  demanda_ruta <- c() # Vector que almacenará la demanda acumulada por cada ruta
  
  # Ruta inicial
  i <- no_asignados[1] #Se toma el primer cliente no asignado (puede escoger el que quiera o escogerlo aleatorio)
  rutas[[1]] <- c(deposito, i, deposito) #Se crea la primera ruta: depósito -> cliente -> depósito
  demanda_ruta <- c(demandas[i]) #Se inicializa la demanda de esa ruta
  no_asignados <- setdiff(no_asignados, i) #Se elimina ese cliente del conjunto de no asignados
  
  while(length(no_asignados) > 0){ #Mientras existan clientes sin asignar
    
    mejor_delta <- Inf #Inicializar el mejor costo de inserción (delta) como infinito
    mejor_cliente <- NULL #Variables para guardar la mejor decisión encontrada, mejor cliente
    mejor_ruta <- NULL #Variables para guardar la mejor decisión encontrada, mejor ruta
    mejor_pos <- NULL #Variables para guardar la mejor decisión encontrada, mejor posición en la ruta
    
    #Evaluar cada cliente no asignado
    for(k in no_asignados){ 
      
      for(r in 1:length(rutas)){ #Evaluar cada ruta existente
        
        ruta <- rutas[[r]] #Extraer la ruta actual
        
        if(demanda_ruta[r] + demandas[k] <= Q){ #Verificar si el cliente cabe en la ruta (restricción de capacidad)
          
          for(pos in 1:(length(ruta)-1)){  #Evaluar todas las posiciones posibles de inserción en la ruta
            
            i <- ruta[pos] #Nodo anterior en la ruta
            j <- ruta[pos+1] #Nodo siguiente en la ruta
            
            delta <- D[i,k] + D[k,j] - D[i,j] #Calcular el costo incremental de insertar k entre i y j
            
            if(delta < mejor_delta){ #Si esta inserción es mejor (menor costo)
              
              #Actualizar la mejor opción encontrada
              mejor_delta <- delta
              mejor_cliente <- k
              mejor_ruta <- r
              mejor_pos <- pos
              
            }
          }
        }
      }
    }
    
    #Si se encontró una inserción factible
    if(!is.null(mejor_cliente)){
      #Extraer la ruta donde se insertará
      ruta <- rutas[[mejor_ruta]]
      #Actualizar la ruta en la lista
      ruta <- append(ruta, mejor_cliente, after = mejor_pos)
      #Actualizar la ruta en la lista
      rutas[[mejor_ruta]] <- ruta
      #Actualizar la demanda acumulada de la ruta
      demanda_ruta[mejor_ruta] <- demanda_ruta[mejor_ruta] + demandas[mejor_cliente]
      #Eliminar el cliente de los no asignados
      no_asignados <- setdiff(no_asignados, mejor_cliente)
      
    } else { #Si no se pudo insertar en ninguna ruta existente
      
      if(length(rutas) < K){ #Verificar si aún hay vehículos disponibles
        
        nuevo <- no_asignados[1] #Tomar un cliente no asignado
        
        rutas[[length(rutas)+1]] <- c(deposito, nuevo, deposito) #Crear una nueva ruta con ese cliente
        
        demanda_ruta <- c(demanda_ruta, demandas[nuevo]) #Registrar su demanda
        
        no_asignados <- setdiff(no_asignados, nuevo) #Eliminarlo de no asignados
        
      } else { #Si no hay más vehículos, el problema no tiene solución factible
        
        stop("No es posible asignar todos los clientes con la flota disponible")
        
      }
    }
  }
  
  return(rutas)
}

Ejecutar algoritmo

#Se ejecuta la función y se obtienen las rutas resultantes.
rutas <- cheapest_insertion(D, demandas, Q, K, deposito)
rutas

Cantidad de rutas

length(rutas)

Cargas de cada ruta

sapply(rutas, function(r) sum(demandas[r]))

Costo de las rutas

distancia_ruta <- function(ruta, D){
  
  distancia <- 0
  
  for(i in 1:(length(ruta)-1)){
    
    distancia <- distancia + D[ruta[i], ruta[i+1]]
    
  }
  
  return(distancia)
}

distancias <- sapply(rutas, distancia_ruta, D)
distancia_total <- sum(distancias)
distancia_total

Resumen de resultados

data.frame(
  Ruta = 1:length(rutas),
  Distancia_km = sapply(rutas, distancia_ruta, D),
  Carga = sapply(rutas, function(r) sum(demandas[r]))
)

Mapa

#Datos
library(readxl)
library(dplyr)
library(leaflet)

#Se lee el archivo de coordenadas geográficas de los municipios.
coords_df <- read_excel("C:/Users/000322041/OneDrive - UPB/UPB/Asignaturas_2026_1/Metodos_cuantitativos_logistica_27043/Clases/3_Heuristica_Insercion_VRP/coordenadas.xlsx")
coords_df$Municipio <- toupper(trimws(coords_df$Municipio)) #Estandarizan los nombres de los municipios
municipios <- toupper(trimws(municipios)) #Estandarizan los nombres del vector municipios


#Preparación de rutas 

#lapply aplica una función a cada ruta para convertirla en coordenadas.
rutas_coords <- lapply(rutas, function(ruta){ #contiene las rutas generadas por Clarke & Wright.
  
  ruta_municipios <- toupper(ruta) #Convierte los nombres de los municipios de la ruta a mayúsculas para garantizar coincidencia con coords_df.
  
  coords_df %>%
    slice(match(ruta_municipios, Municipio)) #Encuentra la posición de cada municipio de la ruta dentro del dataframe de coordenadas. Slece Selecciona esas filas en ese mismo orden.
  
})

#Paleta de colores
colores <- colorFactor(rainbow(length(rutas)), domain = 1:length(rutas)) #Se genera una función de colores para diferenciar rutas. rainbow() crea colores distintos, length(rutas) define cuántos colores se necesitan. 

#Crear mapa: Se crea un mapa interactivo vacío.
mapa <- leaflet() %>%
  addTiles()

#Dibujar cada ruta
for(i in 1:length(rutas_coords)){ #Se recorre cada ruta generada por el algoritmo.
  
  ruta_df <- rutas_coords[[i]] #Se obtiene el dataframe de coordenadas de la ruta i.
  
  mapa <- mapa %>%
    addPolylines(
      lng = ruta_df$longitud,
      lat = ruta_df$latitud,
      color = colores(i),
      weight = 4
    ) %>%
    
    addCircleMarkers(
      lng = ruta_df$longitud,
      lat = ruta_df$latitud,
      radius = 5,
      color = colores(i),
      popup = paste("Vehículo", i, "<br>", ruta_df$Municipio)
    )
}

mapa

Función Relocation inter - ruta

# Relocation inter-ruta (1-0 move)
relocation_interruta <- function(rutas, D, demandas, Q, max_iter = 100){
  
  rutas <- lapply(rutas, function(r) r)  # copia
  mejoro <- TRUE
  iter <- 0
  
  # función auxiliar: carga
  carga_ruta <- function(r){
    sum(demandas[r])
  }
  
  while(mejoro && iter < max_iter){
    
    iter <- iter + 1
    mejoro <- FALSE
    
    for(r1 in 1:length(rutas)){
      
      for(r2 in 1:length(rutas)){
        
        if(r1 == r2) next
        
        ruta1 <- rutas[[r1]]
        ruta2 <- rutas[[r2]]
        
        # recorrer clientes de ruta1 (sin depósito)
        for(i in 2:(length(ruta1)-1)){
          
          cliente <- ruta1[i]
          
          # nueva ruta1 sin cliente
          nueva1 <- ruta1[-i]
          
          # verificar capacidad en ruta2
          nueva_carga2 <- carga_ruta(ruta2) + demandas[cliente]
          if(nueva_carga2 > Q) next
          
          # buscar mejor posición en ruta2
          mejor_delta <- Inf
          mejor_pos <- NULL
          
          for(pos in 1:(length(ruta2)-1)){
            
            a <- ruta2[pos]
            b <- ruta2[pos+1]
            
            delta <- D[a, cliente] + D[cliente, b] - D[a, b]
            
            if(delta < mejor_delta){
              mejor_delta <- delta
              mejor_pos <- pos
            }
          }
          
          # construir nueva ruta2
          nueva2 <- append(ruta2, cliente, after = mejor_pos)
          
          # evaluar mejora total
          costo_actual <- distancia_ruta(ruta1, D) + distancia_ruta(ruta2, D)
          costo_nuevo  <- distancia_ruta(nueva1, D) + distancia_ruta(nueva2, D)
          
          if(costo_nuevo < costo_actual - 1e-6){
            
            rutas[[r1]] <- nueva1
            rutas[[r2]] <- nueva2
            
            mejoro <- TRUE
            break
          }
        }
        
        if(mejoro) break
      }
      
      if(mejoro) break
    }
  }
  
  return(rutas)
}

Ejecución

rutas_reloc <- relocation_interruta(rutas, D, demandas, Q)

Evaluación

distancias_reloc <- sapply(rutas_reloc, distancia_ruta, D)

data.frame(
  Ruta = 1:length(rutas),
  Original = sapply(rutas, distancia_ruta, D),
  Relocation = distancias_reloc,
  Mejora = sapply(rutas, distancia_ruta, D) - distancias_reloc
)

sum(sapply(rutas, distancia_ruta, D))
sum(distancias_reloc)

Cantidad de rutas

length(rutas_reloc)

Cargas de cada ruta

sapply(rutas_reloc, function(r) sum(demandas[r]))