El problema de dispersión de instalaciones modelos de programación de enteros mixtos
Problema de dispersión de instalaciones
Un tutorial completo de Python sobre los modelos de dispersión p y maxisum
En algunos problemas de ubicación de instalaciones, es necesario posicionar las instalaciones de manera que la influencia de una no opaque o afecte negativamente a las demás. Ya sea por motivos de mitigación de riesgos, evitar la competencia u otras consideraciones, dominar las soluciones a estos desafíos es una habilidad crucial en el conjunto de herramientas de investigación de operaciones.
Kuby (1987) propuso dos formulaciones distintas de programación entera mixta (MIP, por sus siglas en inglés) para estos problemas: el problema de dispersión p y el problema de maxisum. En este artículo, ambos serán implementados utilizando la biblioteca de Python Pyomo y el solucionador HiGHS.
Además de los dos modelos, se presentarán algunos recursos de modelado útiles. Primero, una estrategia para linealizar el producto de variables binarias en el contexto de MIP, aunque considerando los objetivos de maximización, no es necesario que sea explícito en este problema. Segundo, una formulación MIP de máx-mín, que pretende maximizar algo menor que cualquier parámetro de un conjunto de elementos si el elemento es seleccionado. Por último, una estrategia para resolver múltiples objetivos con una jerarquía de prioridades bien definida que combina elementos de los dos modelos.
Aquellos interesados que aún no estén familiarizados con la optimización numérica pueden encontrar útil echar un vistazo a mis historias anteriores sobre Programación Lineal y el método Branch & Bound antes de continuar con esta.
Como de costumbre, puedes encontrar el código completo en este repositorio de git.
- Cómo las soluciones GenAI revolucionan la automatización empresaria...
- Novedades en la Cumbre de IA y Negocios e Innovación
- El Congreso llevará a cabo nuevas audiencias sobre la IA
¿Cuál de estas ubicaciones elegirías para colocar 5 instalaciones?
El producto de variables binarias
Cuando se definen los elementos básicos de este problema, se pueden utilizar variables binarias que asumen el valor de 1 si se selecciona una ubicación y 0 en caso contrario. Denotemos estas variables como xᵢ. Supongamos que las distancias (o alguna otra métrica de dispersión) entre dos ubicaciones se calculan de antemano y denotemos estas como dᵢⱼ. ¿Cómo podríamos calcular la dispersión de cualquier par de instalaciones seleccionadas?
De alguna manera intuitiva en esta situación, se puede formular utilizando el producto de las variables binarias xᵢ y xⱼ para calcular la disimilitud obtenida cuando ambas se incluyen en la solución. Sin embargo, esto llevaría a una formulación cuadrática, por lo que no se podrían aplicar solucionadores lineales. Afortunadamente, existe una forma de formular el producto de variables binarias en MIP utilizando algunas restricciones.
Consideremos un grafo dirigido G(V, A) y zᵢⱼ una nueva variable binaria que indica que ambos nodos i y j están seleccionados. Si i o j no están seleccionados, zᵢⱼ debe ser 0. Esto lleva al primer grupo de restricciones:
Hasta ahora, incluso si tanto i como j están seleccionados, zᵢⱼ puede asumir el valor de 0. Por lo tanto, debemos incluir una restricción adicional para que cuando i y j estén seleccionados, zᵢⱼ se convierta en 1.
Cuando se maximiza algo proporcional a zᵢⱼ, como en el problema de maxisum, el segundo grupo de restricciones presentado no es necesario, ya que no tendría sentido no calcular una ganancia proporcional a zᵢⱼ si es factible. Sin embargo, esto puede ser útil al formular otros modelos de MIP.
En la siguiente sección, apliquemos estos conceptos al problema del máximo.
El modelo máximo
El problema máximo discreto debe definir la ubicación de p instalaciones entre N nodos discretos para maximizar la suma de distancias (o el promedio de distancias) calculadas entre todas las parejas de instalaciones seleccionadas. Por lo tanto, consideremos que las instalaciones son nodos desplazados en un grafo dirigido G(V, A). El peso de cada arco de i a j es la métrica de distancia (dispersión) dᵢⱼ conocida de antemano. Además, consideremos variables binarias xᵢ para indicar si se selecciona una ubicación i y zᵢⱼ para indicar si tanto i como j son seleccionados.
El objetivo se puede expresar de la siguiente manera:
Además de las restricciones para linealizar el producto de variables binarias presentadas en la sección anterior, es necesario incluir una restricción para garantizar que se seleccionen p ubicaciones.
Por lo tanto, podemos establecer las restricciones del problema de la siguiente manera:
Vamos a poner eso en código usando Python. Para hacerlo, debemos comenzar importando las bibliotecas que se utilizarán. La biblioteca numpy se utilizará para cálculos de álgebra lineal y para crear puntos de coordenadas aleatorias; las funciones squareform y pdist de scipy se utilizarán para calcular distancias dada una matriz de coordenadas; matplotlib será nuestra herramienta de visualización; y pyomo será el lenguaje de modelado algebraico (AML) (con el solucionador HiGHS).
import numpy as np
from scipy.spatial.distance import squareform, pdist
import matplotlib.pyplot as plt
import pyomo.environ as pyo
from pyomo.contrib.appsi.solvers.highs import HighsDebemos crear N puntos y definir cuántos de ellos deben seleccionarse como ubicaciones de las instalaciones. Al fijar la semilla aleatoria (12) en cada ejecución del código se deben obtener los puntos representados en la introducción.
# Fijar semillas aleatoriasnp.random.seed(12)# Crear puntos aleatoriosN = 25p = 5coordinates = np.random.random((N, 2))# Calcular distancias entre paresdistances = squareform(pdist(coordinates))Ahora tenemos los elementos necesarios para comenzar nuestro modelo pyomo.
Hay dos enfoques para modelar un problema en pyomo: modelos abstractos y concretos. En el primer enfoque, se definen las expresiones algebraicas del problema antes de suministrar algunos valores de datos, mientras que en el segundo, la instancia del modelo se crea inmediatamente a medida que se definen sus elementos. Puede encontrar más información sobre estos enfoques en la documentación de la biblioteca o en el libro de Bynum et al. (2021). A lo largo de este artículo, adoptaremos la formulación del modelo concreto.
model = pyo.ConcreteModel()En este modelo, tenemos dos conjuntos: nodos y arcos. Como estamos considerando un grafo dirigido completo, existen arcos entre cada par de dos nodos excepto desde el nodo hacia sí mismo.
# Conjuntos de nodos y arcosmodel.V = pyo.Set(initialize=range(N))model.A = pyo.Set(initialize=[(i, j) for i in model.V for j in model.V if i != j])Los parámetros proporcionados de antemano son el número de nodos que deben seleccionarse y las distancias entre pares de nodos.
# Parámetrosmodel.d = pyo.Param(model.A, initialize={(i, j): distances[i, j] for (i, j) in model.A})model.p = pyo.Param(initialize=p)Luego, incluyamos las variables de decisión.
# Variables de decisiónmodel.x = pyo.Var(model.V, within=pyo.Binary)model.z = pyo.Var(model.A, within=pyo.Binary)Y las restricciones…
# Se seleccionan p nodosdef p_selection(model): return sum(model.x[:]) == model.p# Si el nodo de inicio no está seleccionado, el arco es 0def dispersion_c1(model, i, j): return model.z[i, j] <= model.x[i]# Si el nodo final no está seleccionado, el arco es 0def dispersion_c2(model, i, j): return model.z[i, j] <= model.x[j]# Incluir restricciones en el modelomodel.p_selection = pyo.Constraint(rule=p_selection)model.dispersion_c1 = pyo.Constraint(model.A, rule=dispersion_c1)model.dispersion_c2 = pyo.Constraint(model.A, rule=dispersion_c2)Finalmente, debemos crear la función objetivo.
def disp_obj(model): return sum(model.z[i, j] * model.d[i, j] for (i, j) in model.A)model.obj = pyo.Objective(rule=disp_obj, sense=pyo.maximize)Ahora el modelo de maximización está listo para ser resuelto. Luego debemos instanciar un solver compatible con pyomo para manejarlo. El solver Highs está disponible en Pyomo (verificar las importaciones) desde la versión 6.4.3 si el paquete highspy está instalado. Así que asegúrese de ejecutar pip install highspy.
solver = Highs()solver.solve(model)Y, después de aproximadamente 120s de tiempo de cómputo, tenemos los siguientes resultados:
Obsérvese que, aunque la dispersión total se maximiza, dos instalaciones en la esquina superior izquierda siguen estando bastante cerca una de la otra, lo cual puede ser indeseable. Por lo tanto, como alternativa a la formulación de maximización, tenemos el modelo de p-dispersión en el cual maximizamos la distancia mínima entre cualquier par de nodos seleccionados.
El modelo de p-dispersión
En el modelo de p-dispersión, necesitamos una variable de decisión adicional para calcular la distancia mínima entre todos los pares de nodos seleccionados, que es nuestro objetivo a maximizar. Denotemos esta variable como D. Debemos crear una restricción de M grande que asegure que si tanto i como j están seleccionados, D sea menor o igual a dᵢⱼ; de lo contrario, debemos asegurar que D sea infinito. Si mantenemos la formulación con zᵢⱼ como el producto de variables binarias, podemos formular esta restricción adicional de la siguiente manera.
En esta formulación, M es un parámetro fijo arbitrariamente grande utilizado para formular una regla disjuntiva. Una buena elección de M debe ser lo suficientemente grande como para garantizar la factibilidad de la restricción (i, j) para cualquier valor de D cuando zᵢⱼ es cero, pero lo más pequeño posible para que la relajación lineal del modelo sea similar a la versión entera, lo que facilita la convergencia. En este problema, el valor más grande de dᵢⱼ puede ser una alternativa interesante.
Aunque esta formulación funcionaría bien para este modelo, se puede utilizar un enfoque más conciso en el cual las variables zᵢⱼ ni siquiera se incluyen en el modelo.
En esta formulación, tanto xᵢ como xⱼ iguales a cero son una condición suficiente para garantizar que la desigualdad sea válida para cualquier valor de D. Y el objetivo se convierte simplemente en maximizar D.
En el siguiente código Python, consideremos que tenemos un nuevo modelo con los mismos conjuntos y parámetros que el anterior, así como el grupo de variables de decisión x.
# Restricción de máximo-mínimo def maxmin_rule(model, i, j): return model.D <= model.d[i, j] + model.M * (1 - model.x[i]) + model.M * (1 - model.x[j])# Nuevo parámetro Mmodel.M = max(model.d[:, :])# Nueva variablemodel.D = pyo.Var(within=pyo.NonNegativeReals)# Nueva restricciónmodel.maxmin_rule = pyo.Constraint(model.A, rule=maxmin_rule)# Objetivomodel.obj = pyo.Objective(expr=model.D, sense=pyo.maximize)Y después de llamar al solver, se tardó menos de 1.2s en obtener los siguientes resultados.
Lo cual parece genial ya que las ubicaciones están bien distribuidas en el espacio.
¿Hay alguna manera de mejorar esta distribución?
Un enfoque multicriterio
Recuerda que la función objetivo del modelo de p-dispersión depende solo de la distancia mínima entre cualquier par de nodos seleccionados. Por lo tanto, se pueden obtener varias soluciones utilizando combinaciones de los dos puntos que definen esta distancia y otros con distancias entre sí mayores o iguales al objetivo en sí. ¿Podemos refinar nuestra solución seleccionando la mejor de estas alternativas? Eso conduce al “enfoque multicriterio” de dos pasos propuesto por Kuby (1987).
En el primer paso, se resuelve el modelo de p-dispersión de manera óptima y se almacena el objetivo actual en un parámetro d_opt. Luego, se resuelve un modelo de maxisum con una restricción adicional que asegura que D ≥ d_opt para obtener, entre las soluciones óptimas del modelo de p-dispersión, aquella con la máxima dispersión total.
Para hacer eso en Python, considera que tienes un modelo de p-dispersión instanciado con variables de decisión y restricciones también del modelo de maxisum.
# D debe ser óptimodef composed_constr(model): return model.D >= model.d_opt# Resolver p-dispersionsolver.solve(model)# Nuevo parámetromodel.d_opt = pyo.Param(initialize=model.obj())# Desactivar objetivo anteriormodel.obj.deactivate()# Más variablesmodel.z = pyo.Var(model.A, within=pyo.Binary)# La solución no empeorará el D actualmodel.composed_constr = pyo.Constraint(rule=composed_constr)# Nuevo objetivomodel.obj_disp = pyo.Objective(rule=disp_obj, sense=pyo.maximize)solver.solve(model)Y esto es sorprendentemente rápido porque cuando el solucionador ingresa al segundo objetivo, el espacio de alternativas factibles es sustancialmente pequeño. En menos de un (segundo) adicional, se obtienen los siguientes resultados.
Lectura adicional
Cuando los clientes no están distribuidos de manera uniforme, las instalaciones tienen una capacidad limitada o quizás no se conoce de antemano el número adecuado de instalaciones, probablemente te enfrentes a un problema de Ubicación de Instalaciones diferente. Puedes encontrar la formulación de Balinski (1965) implementada en Python utilizando PuLP en esta increíble historia de Nicolo Cosimo Albanese.
Optimización: Problema de Ubicación de Instalaciones Capacitado en Python
Encuentra el número óptimo y la ubicación de almacenes para reducir costos y satisfacer la demanda
towardsdatascience.com
Conclusiones
En este artículo, se implementaron dos modelos de programación entera mixta para el problema de dispersión de instalaciones en Python utilizando Pyomo. Como se verificó previamente en la literatura, el modelo de maxisum puede conducir a una distribución desigual de elementos en el espacio, mientras que el modelo de p-dispersión produce soluciones con ubicaciones bien distribuidas en el espacio y uniformemente distribuidas. Un objetivo de maxisum se puede utilizar para refinar soluciones del modelo de p-dispersión seleccionando aquella con la mayor dispersión total del conjunto de soluciones óptimas. El código completo utilizado en estos ejemplos está disponible para su uso posterior.
Referencias
Balinski, M. L. 1965. Integer programming: methods, uses, computations. Management science, 12(3), 253–313.
Bynum, M. L., Hackebeil, G. A., Hart, W. E., Laird, C. D., Nicholson, B. L., Siirola, J. D., … & Woodruff, D. L. 2021. Pyomo-optimization modeling in python (Vol. 67). Berlin/Heidelberg, Germany: Springer.
Kuby, M. J. 1987. Programming models for facility dispersion: The p‐dispersion and maxisum dispersion problems. Geographical Analysis, 19(4), 315–329.
