Trabajo Fin de Máster

Trabajo Fin de Máster Gabriel González Domínguez Trabajo Fin de Máster Gabriel González Domínguez Tutorizada por Prof. Justo Puerto Albandoz Índice general Abstract 1 1. Introducción 3 2. Herramientas 5 2.1. Algoritmos de ramificación y corte . . . . . . . . . . . . . . . . . . . 5 2.1.1. Heurística primaria . . . . . . . . . . . . . . . . . . . . . . . . 8 2.1.2. Preprocesamiento . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.1.3. Familia de planos de corte . . . . . . . . . . . . . . . . . . . . 8 2.1.4. Cuándo agregar planos de corte . . . . . . . . . . . . . . . . . 9 2.1.5. Elevación de cortes . . . . . . . . . . . . . . . . . . . . . . . . 10 2.1.6. Detalles de la implementación . . . . . . . . . . . . . . . . . . 11 2.2. Relajación lagrangiana . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.2.1. Problema de relajación. Ramificación y poda . . . . . . . . . . 12 2.2.2. Técnicas de relajación lagrangiana . . . . . . . . . . . . . . . 15 2.2.3. Relajación lagrangiana y restricciones de desigualdad . . . . . 17 2.2.4. Resolver el problema dual lagrangiano . . . . . . . . . . . . . 18 2.2.5. Ejemplo. Rutas más cortas restringidas . . . . . . . . . . . . . 20 ii trabajo fin de máster 2.2.6. Técnica de optimización subgradiante . . . . . . . . . . . . . 22 2.2.7. Relajación lagrangiana y programación lineal . . . . . . . . . 25 2.3. Generación de columnas . . . . . . . . . . . . . . . . . . . . . . . . . 28 2.3.1. Introducción ........................... 28 2.3.2. Reformulación de Dantzig-Wolfe de un problema entero . . . 30 2.3.3. Resolución del problema maestro . . . . . . . . . . . . . . . . 32 2.3.4. Eficacia del problema maestro lineal . . . . . . . . . . . . . . 34 2.3.5. Resolución del problema del viajante simétrico mediante generación de columnas: un ejemplo ilustrativo . . . . . . . . . 35 3. Problema de flujo de mínimo coste: algoritmo simplex en redes 39 3.1. Introducción................................ 39 3.2. Problemas de flujo con coste mínimo . . . . . . . . . . . . . . . . . . 40 3.3. Soluciones de árbol de expansión y soluciones sin ciclos . . . . . . . . 47 3.4. Manteniendo una estructura de árbol de expansión . . . . . . . . . . 52 3.5. Cálculo de potencial de nodo y flujo en nodo . . . . . . . . . . . . . . 53 3.5.1. Calcular los flujos de arco . . . . . . . . . . . . . . . . . . . . 55 3.6. Algoritmo simplex en redes . . . . . . . . . . . . . . . . . . . . . . . . 57 3.6.1. Obtención de una estructura de árbol de expansión inicial . . 58 3.6.2. Condiciones de optimalidad . . . . . . . . . . . . . . . . . . . 58 3.6.3. Arcos que abandonan el árbol . . . . . . . . . . . . . . . . . . 60 3.6.4. Actualización del árbol . . . . . . . . . . . . . . . . . . . . . . 62 3.6.5. Terminación............................ 63 3.7. Árbol de expansión altamente factible . . . . . . . . . . . . . . . . . . 65 ÍNDICE GENERAL iii 3.7.1. Algoritmo simplex en redes para el problema del camino más corto................................ 70 3.7.2. Algoritmo simplex en redes para el problema de flujo máximo 74 3.7.3. Análisis de sensibilidad . . . . . . . . . . . . . . . . . . . . . . 79 3.7.4. Unimodularidad.......................... 82 4. Problema de flujo multiservicio 85 4.1. Introducción................................ 85 4.1.1. Supuestos............................. 86 4.2. Condiciones de optimalidad . . . . . . . . . . . . . . . . . . . . . . . 87 4.3. Relajación Lagrangiana . . . . . . . . . . . . . . . . . . . . . . . . . . 89 4.4. Generación de columnas . . . . . . . . . . . . . . . . . . . . . . . . . 93 4.4.1. Reformulación de flujos en ruta . . . . . . . . . . . . . . . . . 93 4.4.2. Condiciones de optimalidad . . . . . . . . . . . . . . . . . . . 98 4.5. Procedimiento de solución de generación de columnas . . . . . . . . 99 4.5.1. Determinación de cotas inferiores . . . . . . . . . . . . . . . . 102 4.6. Descomposición de Dantzig-Wolfe . . . . . . . . . . . . . . . . . . . . 103 4.7. Descomposición de la directiva de recursos . . . . . . . . . . . . . . . 105 4.8. Resolver los modelos de directiva de recursos . . . . . . . . . . . . . . 107 4.9. Particiónbase ............................... 110 5. Métodos basados en los flujos en redes para el problema del óptimo enrutamiento de tuberías en diseños navales 117 5.1. Introduction ................................ 117 5.2. El problema de enrutamiento en el contexto del diseño naval. . . . . . 118 iv trabajo fin de máster 5.3. Un modelo basado en el flujo de red multiservicio (Faltan referencias ydibujos).................................. 121 5.3.1. Un método de solución exacta para el (PTTB) . . . . . . . . . 126 5.3.2. Heurística............................. 128 5.4. Experiencia computacional (provisional, habría que cambiarlo) . . . . 129 5.4.1. Instancias aleatorias . . . . . . . . . . . . . . . . . . . . . . . 129 5.4.2. Casodeestudio.......................... 131 6. Conclusión 137 Abstract It is reviewes the theory involved in multicommodity-flow problems and it is analyzed several aspects concerning the optimal routing of pipelines in naval design motivated by a recent collaboration with a leading Naval Engineering company. We start by considerering branch-and-cut algorithm which is a important tool to implement some desired properties in naval design as the minimum required distance between the pipelines volumes or the minimum allowed distances between pipe elbows. We also study Lagrangian relaxation and columns generation, both imperative for a deeper comprehesion of multicommodity-flow problems. In addition, it is studied the minimum cost flows as the the most basic of all flow problems and base for solving multicommodity-flow problems. Moreover, we apply all these mathematical tools into multicommodity-flow problems. Finally, it is proposed a general methodology for the automatic routing. We construct a network-shape framework for feasible solutions of the problem and a cost-based structure for the network that incorporates different desired physical characteristics for the routings. Then, a multicommodity-flow based approach is proposed and a heuristic algorithm is designed to solve it that take into account some desired properties of the routes. 8trabajo fin de máster válido para el primer subproblema. Esta desigualdad se puede modificar para que sea válida para el primer subproblema utilizando una técnica de elevación de coeficientes. 2.1.1 Heurística primaria En el problema de ejemplo 2.1, una vez encontrado un plano de corte apropiado fue posible podar el segundo subproblema por acotaciones. En este caso, la solución de la relajación de la programación lineal del primer subproblema es entera, proporcionando una buena solución. En muchos casos, se necesitan muchas iteracciones hasta que la solución de la relajación es entera. Por lo tanto, a menudo es útil tener buenas heurísticas para convertir la solución fraccional de una relajación en una buena solución entera que se puede utilizar para eliminar otros subproblemas. 2.1.2 Preprocesamiento Un componente muy importante de un algoritmo práctico de ramificación y corte es el preprocesamiento para eliminar restricciones innecesarias, determinar una variable fija y simplificar el problema de otras formas. 2.1.3 Familia de planos de corte Algunas familias de planos de corte son los cortes de Chvátal-Gomory o los cortes derivados de problemas de mochila con variables binarias. Cualquier desigualdad en variables binarias se puede representar como una desigualdad de mochila ∑𝑖∈𝑁𝑎𝑖𝑥𝑖≤ 𝑏con todo 𝑎𝑖>0para un subconjunto 𝑁de variables, donde se eliminan variables o se reemplazan una variable 𝑥𝑗por 1 − 𝑥𝑗cuando sea necesario. La estructura facial del politopo de mochila se puede utilizar para derivar desigualdades válidas para el problema. Por ejemplo, si 𝑅⊆𝑁con ∑𝑖∈𝑅𝑎𝑖> 𝑏, entonces ∑𝑖∈𝑅𝑥𝑖≤|𝑅|− 1 es una desigualdad válida. Además, si 𝑅es un conjunto mínimo (eliminar cualquier elemento de 𝑅hace que ∑𝑖∈𝑅𝑎𝑖< 𝑏), entonces la desigualdad define una cara del politopo de mochila correspondiente. Otras desigualdades pueden derivarse del politopo de mochila. Estas desigualdades se extienden a mochilas con variables enteras generales y una variable continua y a problemas binarios con cotas superiores generalizados. Otra familia de desigualdades útiles son las de elevación y proyección odesigualdades 2. herramientas 9 disyuntivas. Dada la región factible para un problema de programación binaria 𝑆∶= {𝑥∶𝐴𝑥 ≤𝑏, 𝑥𝑖∈ {0,1},∀𝑖}, cada variable puede usarse para generar un conjunto de desigualdades disyuntivas. Sea 𝑆0 𝑗∶= {𝑥∶𝐴𝑥 ≤𝑏, 0≤𝑥𝑖≤1,∀𝑖, 𝑥𝑗= 0} y 𝑆1 𝑗∶= {𝑥∶𝐴𝑥 ≤𝑏, 0≤𝑥𝑖≤1,∀𝑖, 𝑥𝑗= 1}. Entonces 𝑆 ⊆ 𝑆0 𝑗∪𝑆1 𝑗, por lo que se pueden generar desigualdades válidas para 𝑆encontrando desigualdades válidas para la envolvente convexa de 𝑆0 𝑗∪𝑆1 𝑗. Estas desigualdades se generan al resolver problemas de programación lineal. Debido al coste, los cortes generalmente sólo se generan en el nodo raíz. No obstante, pueden ser muy eficaces computacionalmente. 2.1.4 Cuándo agregar planos de corte La sobrecarga computacional de buscar planos de corte puede resultar prohibitiva. Por lo tanto, es común no buscar en todos los nodos del árbol. Una alternativa es buscar en cada octavo nodo o en cada nodo a una profundidad de un múltiplo de ocho en el árbol. Generalmente, en cada nodo del árbol de ramificación y acotación se resuelve la relajación de programación lineal, se encuentran planos de corte, estos se añade a la relajación y se repite el proceso. Se suele llegar a un punto en el que el proceso decae, es decir, la solución de una relajación no es mucho mejor que las soluciones de las relajaciones anteriores. Es entonces aconsejable dejar de trabajar en este nodo y esta rama. Esta reducción gradual es más por la falta de conocimiento sobre la estructura poliédrica de la relajación, que por debilidad del método de planos de corte. En algunas implementaciones, se realiza un número fijo de rondas de búsqueda de planos de corte en un nodo, con quizás varias rondas realizadas en el nodo raíz y un número menor de rondas más abajo en el árbol. La variante de corte y ramificación agrega planos de corte sólo en el nodo raíz del árbol. Por lo general, la implementación de dicho método requerirá un gran esfuerzo para generar planos de corte, requiriendo un tiempo mucho mayor que simplemente resolver la relajación en la raíz. Los beneficios de cortar y ramificar incluyen: Los cortes generados son válidos en la raíz, luego también son válidos en todo el árbol. La contabilidad se reduce, ya que las relajaciones son idénticas en todos los nodo. No se pierde tiempo generando planos de corte en otros nodos. 10 trabajo fin de máster Cortar y ramificar es, en general, una técnica excelente para muchos problemas de programación entera, pero carece del poder que tiene ramificar y cortar en el caso de algunos problemas difíciles. 2.1.5 Elevación de cortes Un corte añadido en un nodo del árbol de ramificación y corte puede no ser válido para otro subproblema. En caso de que el corte no se añada a ningún otro nodo se llama corte local. Este corte sólo afectará al subproblema actual y sus descendientes. El inconveniente de este método es el posible requisito de memoria de tener que almacenar una versión diferente del problema para cada nodo del árbol. Para que un corte sea válido en todo el árbol (corte global), es necesario hacerlo más arriba. La elevación se puede considerar como un método para ajustar una restricción. Volviendo al problema de ejemplo una vez más, la restricción 𝑥1+ 2𝑥2≤6es válida si 𝑥2≤2. Para extender esta restricción y que sea válida cuando 𝑥2≥3, considérese la desigualdad 𝑥1+ 2𝑥2+𝛼(𝑥2− 2) ≤6.(2.15) Se desea elegir un 𝛼tan grande como sea posible mientras sea seguro que se trata de una desigualdad válida. Si 𝑥2= 3 entonces 𝑥1≤1, luego la desigualdad es válida para 𝑥2= 3 siempre que 𝛼≤−1. Si 𝑥2= 1, la desigualdad es válida siempre que 𝛼≥−2. Luego, la desigualdad es válida cuando 𝑥2= 0 siempre que 𝛼≥−2.5. La combinación de estas condiciones da que el rango válido es −2 ≤𝛼≤−1. Las dos opciones extremas 𝛼= −1 y𝛼= −2 dan las desigualdades válidas 𝑥1+𝑥2≤4y 𝑥1≤2, respectivamente. Otras opciones válidas para 𝛼dan desigualdades que son combinaciones convexas de estas dos. De esta manera, las desigualdades válidas para un nodo del árbol se pueden extender a desigualdades válidas en todos los nodos. A menudo no es posible eliminar desigualdades porque la cota superior e inferior de los coeficientes entran en conflicto. Por supuesto, si una desigualdad es válida en el nodo raíz, entonces es válida en todo el árbol, por lo que no es necesario elevarla. A continuación se describe el método de cálculo de coeficientes en el caso de problemas binarios. La desigualdad generada en un nodo del árbol generalmente sólo usará las variables que no están fijas en ese nodo. La elevación se puede utilizar para hacer que la desigualdad sea válida en cualquier nodo del árbol. Es necesario aplicar el proceso de elevación para cada variable que se ha fijado en el nodo, examinando el valor 2. herramientas 11 opuesto para esa variable. Por ejemplo en la desigualdad ∑ 𝑗∈𝐽 𝑎𝑗𝑥𝑗≤ℎpara un subconjunto 𝐽 ⊆ {1,…, 𝑛}(2.16) es válida en un nodo donde 𝑥𝑖se ha fijado a cero, la desigualdad elevada toma la forma ∑ 𝑗∈𝐽 𝑎𝑗𝑥𝑗+𝛼𝑖𝑥𝑖≤ℎ(2.17) para unos escalares 𝛼𝑖. Estos escalarares deben maximizarse para que la desigualdad sea lo más fuerte posible. Maximizar los 𝛼𝑖requiere resolver otro problema de programación entera, por lo que puede ser necesario hacer una aproximación. Este proceso debe aplicarse sucesivamente a cada variable que se haya fijado en el nodo. El orden en el que se examinan las variables puede afectar la desigualdad final, y se pueden obtener otras desigualdades válidas elevando más de una variable a la vez. 2.1.6 Detalles de la implementación Los detalles de la implementación incluye, entre otras cuestiones, la selección de nodos, la selección de variables de ramificación y los requisitos de almacenamiento. Normalmente, un algoritmo de ramificación y acotación almacena la solución en un nodo como una lista de índices de variables básicas. Los algoritmos de ramificación y corte pueden requerir más almacenamiento si los cortes se agregan localmente, porque sería necesario poder recrear la relajación actual en cualquier nodo activo con sólo las restricciones apropiadas. Si los cortes se agregan globalmente, basta con almacenar una única representación del problema. Es posible corregir variables utilizando información sobre costes reducidos y el valor de la mejor solución entera factible conocida. Una vez que las variables se han corregido de esta manera, a menudo es posible corregir variables adicionales utilizando implicaciones lógicas. Para aprovechar al máximo la fijación de variables, la reconstrucción del nodo padre se realiza de la siguiente manera. Una vez que se ha seleccionado un nodo padre, no se divide inmediatamente en dos hijos, sino que se resuelve de nuevo utilizando el algoritmo de los planos de corte. Cuando finaliza el procedimiento de planos de corte, se ha reconstruido el vector de costes reducidos óptimo y se utiliza para fijar variables. Muchas implementaciones de ramificación y corte utilizan un conjunto de cortes. Por lo general, se trata de un conjunto de restricciones que se han generado anteriormente y que no se incluyeron en la relajación o que se eliminaron posteriormente porque 12 trabajo fin de máster ya no parecían estar activas. Es fácil verificar estos cortes para ver si hay violaciones y esto generalmente se hace antes de llamar rutinas de separación más complicadas. El conjunto de cortes también permite reconstruir el nodo padre de manera más eficiente, en parte porque disminuyen las dificultades de la reducción. Una forma de permitir la solución de problemas mucho mayores es utilizar cálculo en paralelo. La naturaleza de los algoritmos de ramificación, y corte y ramificación y acotación hace posible que se realicen en paralelo: típicamente, una relajación de programación lineal se resuelve en un núcleo del ordenador. 2.2 Relajación lagrangiana La relajación lagrangiana es un método que aproxima la solución de un problema difícil mediante la solución de un problema más sencillo. Este problema más sencillo se obtiene del problema original al eliminar una o varias de sus restricciones e introducirlas en la función objetivo. La restricción aparece en la función objetivo como penalización y su peso viene dado por el valor del multiplicador de Lagrange 𝜇. En la práctica, este problema relajado suele ser más fácil de resolver que el problema original. Se trata de un método ampliamente usado en optimización discreta y, en particular, en problemas de redes. Todas las demostraciones de los resultados y proposiciones que aparecen en esta sección se encuentran en el libro publicado por Ahuja,Magnanti yOrlin [3] en 1993. 2.2.1 Problema de relajación. Ramificación y poda Considérese el siguiente problema de programación entera: Minimizar 𝑐𝑥 (2.18) sujeto a 𝑥∈𝐹 . (2.19) En esta formulación 𝐹representa el conjunto factible de soluciones de un problema de programación entera. Esto es, el conjunto de soluciones 𝑥= (𝑥1, 𝑥2,…, 𝑥𝐽)del sistema 𝑥=𝑏, (2.20) 𝑥𝑗= 0 ó1para 𝑗= 1,2,…, 𝐽. (2.21) 2. herramientas 13 Conceptualmente este problema se resuelve de forma trivial: se enumera todas las posibles combinaciones de las variables de decisión. Esto es, de todos los posibles vectores de ceros y unos (𝑥1, 𝑥2,…, 𝑥𝐽)que satisfacen 𝑥=𝑏se elige el que devuel- ve un menor valor de la función objetivo 𝑐𝑥. Debido al explosivo crecimiento de la combinatoria, la aplicabilidad de este procedimiento se reduce a problemas pequeños. ¿Cómo se pueden reducir este cálculo? Supóngase 𝐹=𝐹1∪𝐹2. Por ejemplo, se puede obtener 𝐹1añadiendo la restricción 𝑥1= 0 a𝐹y𝐹2añadiendo 𝑥1= 1. Véase que la solución óptima del conjunto factible 𝐹es la mejor de las soluciones óptimas de 𝐹1 y𝐹2. Supóngase que ya se ha encontrado una solución óptima 𝑥 que minimiza {𝑐𝑥 ∶ 𝑥∈𝐹2}y el valor de su función objetivo es 𝑧(𝑥) = 100. El número de soluciones enteras potenciales en 𝐹1es entonces 2𝐽−1, aún serían prohibitivamente grande todas estas posibilidades, salvo que 𝐽sea pequeño. En lugar de intentar resolver el problema entero en la región factible 𝐹1supóngase que se resuelve una versión relajada del problema (por ejemplo relajando restricciones de integralidad). En general, se obtiene una relajación eliminando algunas restricciones del modelo: por ejemplo, reemplazando las restricciones 𝑥𝑗≥0enteras, por la restricción 𝑥𝑗≥0, o eliminando una o más restricciones de la forma 𝛼𝑥 =𝛽. En la relajación lagrangiana no sólo se eliminan restricciones, también cambia la función objetivo del problema. Sólo se requiere que se relaje algunas de las restricciones del problema y que el valor de la función objetivo de la relajación sea una cota inferior del valor de la función objetivo del problema original. Sea 𝑥′una solución óptima a la relajación, y sea 𝑧(𝑥′)el valor de la función objetivo de esta solución. Se consideran cuatro posibilidades: 1. La solución 𝑥′no existe porque el problema relajado no tiene una solución factible. 2. La solución 𝑥′se encuentra en 𝐹1(aunque se hayan relajado algunas de las restricciones). 3. La solución 𝑥′no está en 𝐹1y el valor de su función objetivo 𝑧(𝑥′)satisface la desigualdad 𝑧(𝑥′)≥𝑧(𝑥) = 100. 4. La solución 𝑥′no se encuentra en 𝐹1y el valor de su función objetivo 𝑧(𝑥′) satisface la desigualdad 𝑧(𝑥′)< 𝑧(𝑥) = 100. Estas cuatro alternativas engloban todos los resultados posibles y son mutuamente excluyentes. Debe ocurrir uno de ellas. 14 trabajo fin de máster En los casos 1 a 3, terminan los cálculos: se ha resuelto el problema original sobre el conjunto 𝐹, aunque no se ha resuelto explícitamente ningún problema entero (asumiendo que se ha obtenido la solución sobre el conjunto 𝐹2sin resolver un programa entero). En el caso 1, dado que la relajación del conjunto 𝐹1está vacío, el conjunto 𝐹1también está vacío, por lo que la solución 𝑥resuelve el problema entero original (totalmente). En el caso 2, como se ha encontrado la solución óptima en la relajación (y por lo tanto un superconjunto) del conjunto 𝐹1, y esta solución está en 𝐹1, también se ha encontrado la mejor solución en 𝐹1; por lo tanto, 𝑥o𝑥′es la solución al problema original (la solución con el menor valor de la función objetivo). Téngase en cuenta que en este caso se ha considerado (enumerado) implícitamente todas las soluciones de 𝐹1en el sentido de que se sabe que ninguna solución en este conjunto es mejor que 𝑥. En el caso 3, la solución 𝑥tiene un valor de la función objetivo que es tan bueno como la mejor solución de una relajación de 𝐹1, por lo que tiene un valor de función objetivo que es tan bueno como cualquier solución en 𝐹1. Por tanto, 𝑥resuelve el problema original. Téngase en cuenta que en el caso 3 se ha utilizado la cota de la función objetivo para eliminar las soluciones en el conjunto 𝐹1procedentes de otras consideraciones. En el caso 4, aún no se ha resuelto el problema original. Se puede resolver el problema de minimizar {𝑐𝑥 ∶𝑥∈𝐹1}mediante algún método directo de programación entera o se puede dividir 𝐹1en dos conjuntos 𝐹3y𝐹4. Por ejemplo, se podría obtener 𝐹3de 𝐹imponiendo 𝑥1= 0 y𝑥2= 0 y obtener 𝐹4estableciendo 𝑥1= 0 y𝑥2= 1. Entonces se podría aplicar cualquier relajación o método directo para los problemas definidos en los conjuntos 𝐹3y𝐹4. En un procedimiento general de ramificación y poda, se divide sistemáticamente la región factible 𝐹en subregiones 𝐹1, 𝐹 2, 𝐹 3,…, 𝐹 𝐾. Sea 𝑥la mejor solución factible (respecto al valor de la función objetivo) que se ha obtenido en cálculos anteriores. Supóngase que para cada 𝑘= 1,2,…, 𝐾,𝐹𝑘está vacío o 𝑥𝑘es una solución de una relajación del conjunto 𝐹𝑘y𝑐𝑥 ≤𝑐𝑥𝑘. Entonces ningún punto en ninguna de las regiones 𝐹1, 𝐹 2, 𝐹 3,…, 𝐹 𝐾puede mejorar el valor de la función objetivo que 𝑥, por lo que 𝑥resuelve el problema de optimización original. Sin embargo, si 𝑐𝑥 > 𝑐𝑥𝑘para alguna región 𝐹𝑘, se necesita subdividir esta región ”ramificando” en algunas de las variables (es decir, dividiendo una subregión en dos estableciendo 𝑥𝑗= 0 o𝑥𝑗= 1 para alguna variable 𝑗y definir así dos nuevas subregiones). Siempre que se satisfaga 𝑐𝑥 ≤𝑐𝑥𝑘para todas las subregiones (o se sabe que están vacías), se ha resuelto el problema original. La finalidad de usar el método de ramificación y poda es encontrar una solución ópti- 2. herramientas 15 ma resolviendo solo una pequeña cantidad de relajaciones. Para hacerlo, se necesitaría obtener buenas soluciones rápidamente y obtener buenas relajaciones de modo que el valor de la función objetivo 𝑧(𝑥𝑘)de la solución 𝑥𝑘a la relajación del conjunto 𝐹𝑘 esté cerca de la solución óptima de 𝐹𝑘. En la práctica, al implementar el procedimiento de ramificación y poda, es necesario tomar muchas decisiones con respecto al orden para elegir las subregiones, las variables a las que se ramificará para cada subregión y los mecanismos (por ejemplo, procedimientos heurísticos) que se podrían utilizar para encontrar "buenas" soluciones viables. También es necesario desarrollar buenas relajaciones que permitan obtener cotas inferiores efectivas (ajustadas): si las cotas inferiores son débiles, los casos 2 y 3 rara vez ocurrirán y el procedimiento de ramificación y poda degenerará en una enumeración total. Por otro lado, si las cotas son muy ajustadas, las relajaciones permitirán eliminar gran parte de la enumeración y desarrollar procedimientos de solución muy eficaces. 2.2.2 Técnicas de relajación lagrangiana Considérese el siguiente modelo de optimización donde 𝑥es el vector de variables de decisión: 𝑧∗= m ın 𝑐𝑥 (2.22) s. a 𝑥=𝑏, (2.23) 𝑥∈𝑋. (2.24) En este modelo (2.22)-(2.24) tanto las función objetivo tanto como las restricciones son lineales. Las variables de decisión 𝑥están restringidas a pertenecer a un conjunto dado 𝑋que, como se verá, a menudo modela redes de flujo. Por ejemplo, el conjunto de restricciones 𝑋= {𝑥∶𝑥=𝑞, 0≤𝑥≤𝑢}podrían ser todas las soluciones factibles para un problema de flujo de red con un vector de oferta/demanda 𝑞. O bien, el conjunto 𝑋puede contener los vectores de incidencia de todos los árboles de expansión de un grafo determinado. A menos que se indique lo contrario, se asume que el conjunto 𝑋es finito (por ejemplo, para problemas de flujo de red, es el conjunto finito de soluciones de árbol de expansión). El procedimiento de relajación lagrangiana utiliza la idea de relajar las restricciones lineales explícitas llevándolas a la función objetivo mediante los multiplicadores de 16 trabajo fin de máster Lagrange 𝜇. Así el problema queda como: Minimizar 𝑐𝑥 +𝜇(𝑥−𝑏)(2.25) 𝑠.𝑎 𝑥 ∈𝑋, (2.26) que es una relajación lagrangiana y se refieren a la función 𝐿(𝜇)=m ın{𝑐𝑥 +𝜇(𝑥−𝑏) ∶ 𝑥∈𝑋},(2.27) como la función lagrangiana. Dado que al hacer la relajación lagrangiana se eliminan las restricciones 𝑥−𝑏, la solución del subproblema lagrangiano no necesita ser factible para el problema original (2.22)-(2.24). ¿Es posible obtener alguna información útil sobre el problema original incluso cuando la solución al subproblema de Lagrange no es factible? Lema 2.1 (Principio delimitador de la lagrangiana.).Para cualquier vector 𝜇de los multiplicadores lagrangianos, el valor de la función lagrangiana 𝐿(𝜇)es una cota inferior del valor óptimo de la función objetivo 𝑧∗del problema de optimización original (2.22)-(2.24). Para obtener la mejor cota inferior es necesario resolver el siguiente problema de optimización. 𝐿∗= m ax 𝜇𝐿(𝜇)(2.28) que es el problema dual lagrangiano asociado al problema de optimización original (2.22)-(2.24). Propiedad 2.1 ((Dualidad débil).).El valor óptimo de la función objetivo 𝐿∗del problema dual lagrangiano es siempre una cota inferior del valor óptimo de la función objetivo del problema (2.22)-(2.24) (es decir, 𝐿∗≤𝑧∗). Esto proporciona cotas válidos para comparar los valores de la función objetivo del problema dual lagrangiano y la optimalidad de (2.22)-(2.24) con cualquier elección de los multiplicadores de Lagrange 𝜇y cualquier solución factible 𝑥: 𝐿(𝜇)≤𝐿∗≤𝑧∗≤𝑐𝑥. (2.29) Estas desigualdades garantizan la optimalidad de la solución del problema dual Lagrangiano o del problema original (2.22)-(2.24). 2. herramientas 17 Propiedad 2.2 (Prueba de optimalidad.).(a) Supóngase que 𝜇es un vector de multiplicadores lagrangianos y 𝑥es una solución factible al problema de optimización (2.22)-(2.24) que satisface la condición 𝐿(𝜇) = 𝑐𝑥.𝐿(𝜇)es entonces una solución óptima del problema dual lagrangiano (𝐿∗=𝐿(𝜇)) y 𝑥es una solución óptima del problema de optimización (2.22)-(2.24). (b) Si para alguna elección del vector multiplicador lagrangiano 𝜇, la solución 𝑥∗ de la relajación lagrangiana es factible en el problema de optimización (2.22)- (2.24), entonces 𝑥∗es una solución óptima para el problema de optimización (2.22)-(2.24) y 𝜇es una solución óptima para el problema dual lagrangiano. Obsérvese que por el apartado (b) de la Propiedad 2.2, 𝐿(𝜇) = 𝑐𝑥∗+𝜇(𝑥∗−𝑏) y𝑥∗=𝑏. Por lo tanto, 𝐿(𝜇) = 𝑐𝑥∗y el apartado (a) implica que 𝑥∗resuelve el problema (2.22)-(2.24) y 𝜇resuelve el problema dual lagrangiano. Como se indica en la Propiedad 2.2, el principio de delimitación proporciona una garantía (𝐿(𝜇) = 𝑐𝑥) de que la solución factible 𝑥al problema de optimización (2.22)- (2.24) es una solución óptima. Incluso si 𝐿(𝜇)< 𝑐𝑥, tener la cota inferior permite medir qué tan lejos está una solución de la solución optima. Por ejemplo, si [𝑐𝑥 − 𝐿(𝜇)]∕𝐿(𝜇)≤0.05 se sabe que el valor de la función objetivo de la solución factible 𝑥es como mucho un 5 % de lo óptima que podría ser. Este tipo de cota es muy útil en la práctica: permite evaluar el grado de subóptimalidad de soluciones dadas y permite terminar la búsqueda de una solución óptima cuando se tiene una solución lo suficientemente cerca de la optimalidad. 2.2.3 Relajación lagrangiana y restricciones de desigualdad En el modelo de optimización (2.22)-(2.24), las restricciones: 𝑥=𝑏son todas restricciones de igualdad. En la práctica, a menudo se encuentran modelos, como el problema de la ruta más corta restringida, que se formulan de forma más natural en forma de desigualdad: 𝑥≤𝑏. El problema dual lagrangiano para estos problemas es una ligera variante del anterior: 𝐿∗= m ax 𝜇≥0𝐿(𝜇)(2.30) El único cambio en el problema dual lagrangiano es que los multiplicadores están restringidos a ser no negativos. Sin embargo, existe una diferencia sustancial entre relajar las restricciones de igualdad y las restricciones de desigualdad. Cuando se relajan las restricciones de desigualdad 24 trabajo fin de máster Eligiendo 𝜃𝑘= 1∕𝑘se satisface estas condiciones que aseguran la convergencia del algoritmo en una solución óptima del problema dual. Una variante importante del procedimiento de optimización subgradiente es una adaptación del ”método de Newton” para resolver sistemas de ecuaciones no lineales. Supóngase, como antes, que 𝐿(𝜇𝑘) = 𝑐𝑥𝑘+𝜇𝑘(𝑥𝑘−𝑏); es decir, 𝑥𝑘resuelve el subproblema lagrangiano cuando 𝜇=𝜇𝑘. Supóngase que 𝑥𝑘continúa resolviendo el subproblema lagrangiano a medida que se varía 𝜇; o, dicho de otra manera, se hace una aproximación lineal de 𝑟(𝜇) = 𝑐𝑥𝑘+𝜇(𝑥𝑘−𝑏)a𝐿(𝜇). Supóngase, además, que se conoce el valor óptimo 𝐿∗del problema dual lagrangiano (el cual no se conoce). Entonces uno podría moverse en la dirección del subgradiente hasta que el valor de la aproximación lineal sea exactamente igual a 𝐿∗. En general, se establece una longitud de paso 𝜃𝑘para que 𝑟(𝜇𝑘+1) = 𝑐𝑥𝑘+𝜇𝑘+1(𝑥𝑘−𝑏) = 𝐿∗,(2.44) o dado que, 𝜇𝑘+1 =𝜇𝑘+𝜃𝑘(𝑥𝑘−𝑏), 𝑟(𝜇𝑘+1) = 𝑐𝑥𝑘+ [𝜇𝑘+𝜃𝑘(𝑥𝑘−𝑏)](𝑥𝑘−𝑏) = 𝐿∗.(2.45) Reordenando términos, recordando que 𝐿(𝜇𝑘) = 𝑐𝑥𝑘+𝜇(𝑥𝑘−𝑏), y tomando ‖𝑦‖= (∑𝑗𝑦2 𝑗)1∕2 (norma euclídea del vector 𝑦), se puede resolver para hallar la longitud de paso: 𝜃𝑘=𝐿∗−𝐿(𝜇𝑘) ‖𝑥𝑘−𝑏‖2.(2.46) Dado que no se conoce el valor óptimo de la función objetivo 𝐿∗del problema dual lagrangiano (después de todo, eso es lo que se trata de encontrar), suele utilizarse la siguiente heurística para seleccionar la longitud de paso: 𝜃𝑘=𝜆𝑘[𝑈𝐵 −𝐿(𝜇𝑘)] ‖𝑥𝑘−𝑏‖2.(2.47) En esta expresión, 𝑈𝐵 es una cota superior del valor óptimo de la función objetivo 𝑧∗del problema (2.30) y, por lo tanto, una cota superior de 𝐿∗; y 𝜆𝑘es un escalar elegido (estrictamente) entre 0 y 2. Inicialmente , la cota superior es el valor de la función objetivo de cualquier solución factible conocida del problema (2.30). A medida que avanza el algoritmo, si se genera una solución factible mejor (es decir, de menor 2. herramientas 25 coste) y se utiliza el valor de la función objetivo de esta solución en lugar de la cota superior 𝑈𝐵. En realidad, se eligen los escalares 𝜆𝑘comenzando con 𝜆𝑘= 2 y luego reduciendo 𝜆𝑘en un factor de 2 siempre que el mejor valor de la lagrangiana encontrado hasta ahora no haya aumentado en un número específico de iteraciones. Dado que esta versión del algoritmo no tiene criterios de detención convenientes, generalmente termina después de haber realizado un número específico de iteraciones. Si se aplica la relajación lagrangiana a un problema con restricciones 𝑥≤𝑏expresadas en forma de desigualdad en lugar de las restricciones de igualdad, los multiplicadores de Lagrange 𝜇están restringidos a ser no negativos. La fórmula actualizadora 𝜇𝑘+1 =𝜇𝑘+𝜃𝑘(𝑥𝑘−𝑏)puede hacer que uno o más de las componentes 𝜇𝑖de 𝜇se vuelvan negativas. Para evitar esta posibilidad, se modifica la fórmula actualizadora de la siguiente manera: 𝜇𝑘+1 =[𝜇𝑘+𝜃𝑘(𝑥𝑘−𝑏)]+.(2.48) En esta expresión, la notación [𝑦]+denota la "parte positiva" del vector 𝑦; es decir, el 𝑖-ésimo componente de [𝑦]+es igual al máximo entre 0 e 𝑦𝑖. Dicho de otra manera, si la fórmula actualizadora 𝜇𝑘+1 =𝜇𝑘+𝜃𝑘(𝑥𝑘−𝑏)causase que el 𝑖-ésimo componente de 𝜇𝑖fuera negativo, entonces simplemente se establecería el valor de este componente como cero. Luego, se implementan todos los demás pasos del procedimiento de subgradiante (es decir, la elección del tamaño de paso 𝜃en cada paso y la solución de los subproblemas de Lagrange) exactamente igual que para los problemas con restricciones de igualdad. Para problemas con restricciones de igualdad y desigualdad, se usa una mezcla sencilla de las versiones de igualdad y desigualdad del algoritmo: siempre que la fórmula actualizadora para los multiplicadores de Lagrange cause que cualquier componente 𝜇𝑖de 𝜇correspondiente a una restricción de desigualdad sea negativa, se establece el valor de ese multiplicador en cero. 2.2.7 Relajación lagrangiana y programación lineal Teorema 2.2. Supóngase que se aplica la técnica de relajación lagrangiana a un problema de programación lineal (𝑃′)definido como m ın{𝑐𝑥 ∶𝑥=𝑏, 𝑥≤𝑞, 𝑥 ≥0} relajando las restricciones 𝑥=𝑏. Entonces, el valor óptimo 𝐿∗del problema dual de Lagrange es igual al valor óptimo de la función objetivo de (𝑃′). El Teorema 2.2 anterior muestra que la técnica de relajación lagrangiana proporciona un método alternativo para resolver un problema de programación lineal. En lugar 26 trabajo fin de máster de resolver el problema directamente usando cualquier algoritmo de programación lineal, se puede relajar un subconjunto de las restricciones y resolver el problema dual de Lagrange utilizando la optimización subgradiente y resolver una secuencia de problemas relajados. En algunas situaciones, el problema relajado es fácil de resolver, mientras el problema original no lo es. En estas situaciones, un algoritmo basado en la relajación lagrangiana es un método de solución atractivo. A continuación, supóngase que se aplica la relajación lagrangiana a un problema de optimización discreto (𝑃)definido como m ın{𝑐𝑥 ∶𝑥=𝑏, 𝑥 ∈𝑋}. Sea el conjunto discreto 𝑋= {𝑥∶𝑥≤𝑞, 𝑥 ≥0y entero}para una matriz entera y un vector de enteros 𝑞. En consecuencia, el problema (𝑃)se convierte en: 𝑧∗= m ın{𝑐𝑥 ∶𝑥=𝑏, 𝑥≤𝑞, 𝑥 ≥0y entero}(P) No se incurre en pérdida de generalidad al especificar el conjunto 𝑋de esta manera porque se puede formular casi todos los problemas de optimización discreta de la vida real como problemas de programación entera. Sea (𝐿𝑃 )la relajación de programación lineal del problema (𝑃)y sea 𝑧0el valor óptimo de su función objetivo. Es decir, 𝑧0= m ın{𝑐𝑥 ∶𝑥=𝑏, 𝑥≤𝑞, 𝑥 ≥0} (LP) El conjunto de soluciones factibles de (P) se encuentra dentro del conjunto de soluciones factibles de (LP), luego 𝑧0≤𝑧∗. Por lo tanto, la relajación del problema lineal proporciona una cota inferior válida del valor óptimo de la función objetivo de (P). En la Propiedad 2.1 se ha mostrado anteriormente que el problema dual de Lagrange también proporciona una cota inferior 𝐿∗del valor óptimo de la función objetivo de (P). Ahora se muestra que 𝑧0≤𝐿∗; es decir, la relajación lagrangiana proporciona una cota inferior que es al menos tan buena como la obtenida de la relajación del problema lineal. Se establece este resultado mostrando que el problema dual de Lagrange también resuelve un problema de programación lineal pero el espacio de solución de este problema está contenido dentro del espacio de solución del problema (LP). El problema de programación lineal que resuelve el problema dual de Lagrange utiliza la ”convexificación” del espacio de solución 𝑋= {𝑥∶𝑥≤𝑞, 𝑥 ≥0y entero}. Sea 𝑋= {𝑥1, 𝑥2,…, 𝑥𝐾}un conjunto finito. Una solución 𝑥es una combinación convexa de soluciones 𝑥1, 𝑥2,…, 𝑥𝐾si 𝑥=∑𝐾 𝑘=1 𝜆𝑘𝑥𝑘para pesos no negativos 𝜆1, 𝜆2,…, 𝜆𝐾 que satisfacen la condición ∑𝐾 𝑘=1 𝜆𝑘= 1. Sea (𝑋)la envolvente convexa de 𝑋(es decir, el conjunto de todas las combinaciones convexas de 𝑋). Se tienen las siguientes propiedades de (𝑋). 2. herramientas 27 Propiedad 2.4.1. El conjunto (𝑋)es un poliedro: puede expresarse como un conjunto definido por un número finito de desigualdades lineales. 2. Cada punto extremo de la solución del poliedro (𝑋)está situado en 𝑋, y si se optimiza una función objetivo lineal sobre (𝑋), alguna solución en 𝑋será la solución óptima. 3. El conjunto (𝑋)está contenido en el conjunto de soluciones {𝑥∶𝑥≤𝑞, 𝑥 ≥ 0}. Propiedad 2.5.El valor óptimo 𝐿∗de la función objetivo del problema dual lagrangiano es igual al valor óptimo de la función objetivo del programa lineal m ın{𝑐𝑥 ∶ 𝑥=𝑏, 𝑥 ∈(𝑋)}. La versión convexa del problema (𝑃) m ın{𝑐𝑥 ∶𝑥=𝑏, 𝑥 ∈(𝑋)} es y se le llama (𝐶𝑃 ). El teorema anterior 2.2 muestra que 𝐿∗es igual al valor óptimo de la función objetivo del problema convexo. ¿Cuál es la relación entre el conjunto de soluciones factibles del problema convexo (𝐶𝑃 )y la relajación del problema lineal (𝐿𝑃 )? Obsérvese que de la Propiedad 2.4 se deduce que (𝑋)está contenido en el conjunto {𝑥∶𝑥≤𝑞, 𝑥 ≥0}, el conjunto de soluciones del problema (𝐶𝑃 )dado por {𝑥∶ 𝑥=𝑏, 𝑥 ∈(𝑋)} está contenido en el conjunto de soluciones de (𝐿𝑃 )dado por {𝑥∶𝑥=𝑏, 𝑥≤𝑞, 𝑥 ≥0}. Puesto que optimizar la misma función objetivo en un espacio de solución más pequeño no puede mejorar el valor de la función objetivo, se tiene que 𝑧0≤𝐿∗. Propiedad 2.6.En un problema entero expresado en forma de minimización, la cota inferior obtenida por la relajación lagrangiana es siempre tan grande como la cota obtenida por la relajación del problema lineal; es decir, 𝑧0≤𝐿∗. ¿En qué situaciones la cota de la relajación lagrangiana será igual a la cota de programación lineal? Si el subproblema lagrangiano satisface la propiedad de integralidad, la cota lagrangiana será igual a la cota de programación lineal. Se dice que el subproblema de Lagrange m ın{𝑑𝑥 ∶𝑥≤𝑞, 𝑥 ≥0y entero}satisface la propiedad de integralidad si tiene una solución óptima entera para cada elección de coeficientes de la función objetivo incluso si se relajan las restricciones de integralidad sobre las variables 𝑥. Téngase en cuenta que esta condición supone que los problemas m ın{𝑐𝑥+𝜇(𝑥−𝑏) ∶ 𝑥≤𝑞, 𝑥 ≥0y entero}ym ın{𝑐𝑥+𝜇(𝑥−𝑏) ∶ 𝑥≤𝑞, 𝑥 ≥0} tienen los mismos valores óptimos de la función objetivo para cada elección del multiplicador de Lagrange 𝜇. Por ejemplo, si las restricciones 𝑥≤𝑞son las restricciones de balance de flujo de un problema de flujo de coste mínimo (o cualquiera de sus casos especiales, como el flujo máximo, la ruta más corta y los problemas de asignación), el 28 trabajo fin de máster problema m ın{𝑐𝑥 +𝜇(𝑥−𝑏) ∶ 𝑥≤𝑞, 𝑥 ≥0} siempre tendrá una solución óptima entera y la imposición de restricciones de integralidad a las variables no aumentará el valor óptimo de la función objetivo. Propiedad 2.7.Si el subproblema lagrangiano del problema de optimización (𝑃)satisface la propiedad de integralidad, entonces 𝑧0=𝐿∗. Para problemas que satisfacen la propiedad de integralidad, resolver el problema dual lagrangiano es equivalente a resolver la relajación de programación lineal del problema. En estas situaciones, la técnica de relajación lagrangiana no proporciona una cota mejor que la relajación de programación lineal. No obstante, la técnica de relajación lagrangiana aún podría tener un valor considerable, porque resolver el problema dual lagrangiano podría ser más eficiente que resolver la relajación de programación lineal directamente. Los problemas de optimización de redes quizás conformen el dominio de problemas más útil para explotar este resultado porque el subproblema lagrangiano en estos casos a menudo resulta ser un problema de flujo de coste mínimo o una de sus variantes. En la mayoría de los casos, el valor óptimo de la función objetivo 𝐿∗del problema dual de Lagrange será estrictamente menor que el valor óptimo de la función objetivo 𝑧∗del problema (𝑃). 2.3 Generación de columnas 2.3.1 Introducción Sea un problema de programación entera m ax{𝑐𝑥 ∶𝑥∈𝑋}con una región factible 𝑋que puede escribirse como la intersección de dos o más conjuntos con estructura 𝑋= ∩𝐾 𝑘=0𝑋𝑘para un 𝐾≥1. Tómese el caso en el que las restricciones tienen la forma: 𝐴1𝑥1+𝐴2𝑥2+⋯+𝐴𝐾𝑥𝐾=𝑏(2.49) 𝐷1𝑥1≤𝑑1(2.50) ⋮ 𝐷𝐾𝑥𝐾≤𝑑𝐾(2.51) 𝑥1∈𝑍𝑛1 +,…, 𝑥𝐾∈𝑍𝑛𝐾 +(2.52) 2. herramientas 29 de manera que los conjuntos 𝑋𝑘= {𝑥𝑘∈𝑍𝑛𝑘 +∶𝐷𝑘𝑥𝑘≤𝑑𝑘}son independientes para 𝑘= 1,…, 𝐾, y sólo las restricciones ∑𝐾 𝑘=1 𝐴𝑘𝑥𝑘=𝑏asocian entre sí los diferentes conjuntos de variables. Dada la función objetivo m ax ∑𝐾 𝑘=1 𝑐𝑘𝑥𝑘, dos herramientas que explotan con éxito esta estructura son la generación de cortes, con la cual se trata de validar desigualdades para los subconjuntos 𝑋𝑘,𝑘= 1,…, 𝐾; y la relajación lagrangiana, con la que se dualizan las restricciones ∑𝐾 𝑘=1 𝐴𝑘𝑥𝑘=𝑏para obtener el problema dual: m ın 𝑢𝐿(𝑢),(2.53) donde 𝐿(𝑢)=m ax{ 𝐾 ∑ 𝑘=1 (𝑐𝑘−𝑢𝐴𝑘)𝑥𝑘+𝑢𝑏 ∶𝑥𝑘∈𝑋𝑘para 𝑘= 1,…, 𝐾},(2.54) y el cálculo de 𝐿(𝑢)se divide en 𝐾subproblemas distintos: 𝐿(𝑢) = 𝐾 ∑ 𝑘=1 m ax{(𝑐𝑘−𝑢𝐴𝑘)𝑥𝑘∶𝑥𝑘∈𝑋𝑘} + 𝑢𝑏. (2.55) A continuación, se presenta una tercera herramienta que explota la forma de programación entera anterior. El procedimiento esencialmente resuelve un problema equivalente de la siguiente forma: m ax{ 𝐾 ∑ 𝑘=1 𝛾𝑘𝜆𝑘∶ 𝐾 ∑ 𝑘=1 𝐵𝑘𝜆𝑘=𝛽, 𝜆𝑘≥0entero para 𝑘= 1,…𝐾}(2.56) donde cada matriz 𝐵𝑘tiene un gran número de columnas, una por cada punto factible en 𝑋𝑘, y cada vector 𝜆𝑘contiene las correspondientes variables. Todas las demostraciones de los resultados y proposiciones que aparecen en esta sección se encuentran en el libro publicado por Wolsey [4] en 1998. Véase, por ejemplo, una formulación alternativa de este tipo para el problema de ubicación de instalaciones no capacitadas. Las localizaciones 𝑗= 1,…, 𝑛 se corresponden con los índices 𝑘= 1,…, 𝐾. Si el depósito 𝑗satisface la demanda del conjunto cliente 𝑆, entonces se establece 𝜆𝑗 𝑆= 1 por cada subconjunto no vacío 𝑆 ⊆ 𝑀 de 30 trabajo fin de máster clientes. Esto conduce a la siguiente formulación: m ın ∑ 𝑗∈𝑁∑ 𝑆≠𝜙(∑ 𝑖∈𝑆 𝑐𝑖𝑗 +𝑓𝑗)𝜆𝑗 𝑆(2.57) ∑ 𝑗∈𝑁∑ 𝑆≠𝜙, 𝑖∈𝑆 𝜆𝑗 𝑆= 1 para 𝑖∈𝑀(2.58) ∑ 𝑆≠𝜙 𝜆𝑗 𝑆≤1para 𝑗∈𝑁(2.59) 𝜆𝑗 𝑆∈ {0,1} para 𝜙≠𝑆 ⊆ 𝑀, 𝑗 ∈𝑁. (2.60) El coste 𝜆𝑗 𝑆es el coste de abrir el depósito 𝑗y atender a los clientes en 𝑆desde el depósito 𝑗. Las restricciones (2.58) imponen que se atienda a cada cliente, mientras que las restricciones (2.59) aseguran que al deposito 𝑗se le asigna como mucho un subconjunto de clientes. En la práctica las restricciones (2.60) son tipicamente innecesarias. Luego los problemas que se desean resolver son problemas enteros con un número enorme de variables, donde las columnas a menudo se describen implicitamente como los vectores de incidencia de ciertos subconjuntos de un conjunto que recorre subconjuntos clientes. Estas grandes formulaciones de un problema entero se denominan problemas maestros. A continuación, se considera cómo resolver la relajación de los problemas maestros y comparar esta relajación con la obtenida por la relajación lagrangiana, o por el uso de planos de corte. Finalmente, se considera qué hacer cuando la solución del problema lineal no es entera, y se debe recurrir a la enumeración lo que lleva a un problema entero de generación de columnas o a algoritmos de ramificación y acotación. 2.3.2 Reformulación de Dantzig-Wolfe de un problema entero Considérese el problema entero de la forma: 𝑥= m ax{ 𝐾 ∑ 𝑘=1 𝑐𝑘𝑥𝑘∶ 𝐾 ∑ 𝑘=1 𝐴𝑘𝑥𝑘=𝑏, 𝑥𝑘∈𝑋𝑘para 𝑘= 1,…, 𝐾}(2.61) 2. herramientas 31 donde 𝑋𝐾= {𝑥𝑘∈𝑍𝑛𝑘 +;𝐷𝑘𝑥𝑘≤𝑑𝑘}para 𝑘= 1,…, 𝐾. Asumíendo que cada conjunto 𝑋𝑘contiene un enorme aunque finito conjunto de puntos {𝑥𝑘,𝑡}𝑇𝑘 𝑡=1, se tiene que 𝑋𝑘=⎧ ⎪ ⎨ ⎪ ⎩ 𝑥𝑘∈𝑅𝑛𝑘∶𝑥𝑘=∑𝑇𝑘 𝑡=1 𝜆𝑘,𝑡𝑥𝑥,𝑡, ∑𝑇𝑘 𝑡=1 𝜆𝑘,𝑡 = 1, 𝜆𝑘,𝑡 ∈ {0,1} para 𝑡= 1,…, 𝑇𝑘⎫ ⎪ ⎬ ⎪ ⎭ .(2.62) Sustituyendo ahora por 𝑥𝑘se llega al problema entero maestro equivalente: 𝑧= m ax 𝐾 ∑ 𝑘=1 𝑇𝑘 ∑ 𝑡=1 (𝑐𝑘𝑥𝑘,𝑡)𝜆𝑘,𝑡 (2.63) 𝐾 ∑ 𝑘=1 𝑇𝑘 ∑ 𝑡=1 (𝐴𝑘𝑥𝑘,𝑡)𝜆𝑘,𝑡 =𝑏(2.64) 𝑇𝑘 ∑ 𝑡=1 𝜆𝑘,𝑡 = 1 para 𝑘= 1,…, 𝐾 (2.65) 𝜆𝑘,𝑡 ∈ {0,1} para 𝑡= 1,…, 𝑇𝐾y𝑘= 1,…, 𝐾. (2.66) Retomando el problema de ubicación de instalaciones no capacitadas, supóngase que se parte de la formulación débil: m ın ∑ 𝑖∈𝑀∑ 𝑗∈𝑁 𝑐𝑖𝑗𝑥𝑖𝑗 +∑ 𝑗∈𝑁 𝑓𝑗𝑦𝑗(2.67) ∑ 𝑗∈𝑁 𝑥𝑖𝑗 = 1 para 𝑖∈𝑀(2.68) ∑ 𝑖∈𝑀 𝑥𝑖𝑗 ≤𝑚𝑦𝑗para 𝑗∈𝑁(2.69) 𝑥∈𝐵|𝑀|×|𝑁|, 𝑦 ∈𝐵|𝑁|.(2.70) Donde: 𝑋𝑘= {(𝑥1𝑘,…, 𝑥𝑚𝑘, 𝑦𝑘) ∶ ∑ 𝑖∈𝑀 𝑥𝑖𝑘 ≤𝑚𝑦𝑘, 𝑥𝑖𝑘 ∈𝐵1para 𝑖∈𝑀, 𝑦𝑘∈ {0,1}}. (2.71) 32 trabajo fin de máster Los puntos en 𝑋𝑘son {(𝑥𝑘 𝑆,1)}𝑆⊆𝑀 donde 𝑥𝑘 𝑆es el vector de incidencia de 𝑆 ⊆ 𝑀, y (0,0) con variables asociadas 𝜆𝑘 𝑆,𝜈𝑘que conducen al problema maestro entero: m ın ∑ 𝑗∈𝑁[∑ 𝑆≠𝜙(∑ 𝑖∈𝑆 𝑐𝑖𝑗 +𝑓𝑖)𝜆𝑗 𝑆+𝑓𝑗𝜆𝑗 ∅](2.72) ∑ 𝑗∈𝑁∑ 𝑆≠𝜙, 𝑖∈𝑆 𝜆𝑗 𝑆= 1 para 𝑖∈𝑀(2.73) ∑ 𝑆≠𝜙 𝜆𝑗 𝑆+𝜆𝑗 ∅𝜈𝑗= 1 para 𝑗∈𝑁(2.74) 𝜆𝑗 𝑆∈ {0,1} para 𝑆 ⊆ 𝑀, 𝑗 ∈𝑁, 𝜈𝑗∈ {0,1} para 𝑗∈𝑁. (2.75) Obsérvese que 𝑓𝑗≥0, luego 𝜈𝑗es dominante sobre la variable 𝜆𝑗 ∅y𝜆𝑗 ∅puede eliminarse. Si se coge 𝜈𝑗como variable de holgura en (2.59) entonces las formulaciones (2.57)-(2.60) y (2.72)-(2.75) son idénticas. 2.3.3 Resolución del problema maestro Para resolver la relajación lineal del problema maestro entero se utiliza un algoritmo de generación de columnas llamado problema maestro de programación lineal. 𝑧𝑃 𝑀𝐿 = m ax 𝐾 ∑ 𝑘=1 𝑇𝑘 ∑ 𝑡=1 (𝑐𝑘𝑥𝑘,𝑡)𝜆𝑘,𝑡 (2.76) 𝐾 ∑ 𝑘=1 𝑇𝑘 ∑ 𝑡=1 (𝐴𝑘𝑥𝑘,𝑡)𝜆𝑘,𝑡 =𝑏(2.77) 𝑇𝑘 ∑ 𝑡=1 𝜆𝑘,𝑡 = 1 para 𝑘= 1,…, 𝐾 (2.78) 𝜆𝑘,𝑡 ≥0para 𝑡= 1,…, 𝑇 𝑘, 𝑘 = 1,…, 𝐾 (2.79) donde hay una columna ⎛⎜⎜⎝ 𝑐𝑘𝑥 𝐴𝑘𝑥 𝑒𝑘⎞⎟⎟⎠ para cada 𝑥∈𝑋𝑘. Como variables duales asociadas a las restricciones (2.77) se usará {𝜋𝑖}𝑚 𝑖=1 y para las restricciones (2.78), llamadas restricciones de convexidad, se usarán como variables duales {𝜇𝑘}𝐾 𝑘=1. La idea es resolver el problema lineal con el algoritmo simplex primal. Sin embargo, el paso de tasar y elegir una columna que entre en la base debe modificarse debido al 2. herramientas 33 enorme número de columnas. En lugar de tasar las columnas una a una, se resuelve el problema de encontrar la columna con el mayor coste reducido, que es en sí mismo un conjunto de 𝐾problemas de optimización. Inicialización. Supóngase que se dispone de un subconjunto de columnas (al menos una por cada 𝑘), el cual forma un problema maestro lineal restringido factible: 𝑧𝑃 𝑀𝐿𝑅 = m ax 𝑐  𝜆(2.80)  𝐴 𝜆= 𝑏(2.81)  𝜆≥0(2.82) donde  𝑏=(𝑏 1), el conjunto de columnas disponible genera la submatriz:  𝐴=⎛⎜⎜⎜⎜⎜⎝ 𝐴1𝑥1,1…𝐴1𝑥1,𝑇1𝐴2𝑥2,1…𝐴1𝑥1,𝑇2…𝐴𝐾𝑥𝐾,1…𝐴𝐾𝑥𝐾,𝑇𝐾 1 … 1 1 … 1 ⋱ 1 … 1 ⎞⎟⎟⎟⎟⎟⎠ (2.83) y𝑐, 𝜆son los costes y variables correspondientes. Al resolver (2.80)-(2.82) se obtiene una solución optima primal  𝜆∗y una solución óptima dual (𝜋, 𝜇) ∈ 𝑅𝑚×𝑅𝐾. Factibilidad primal. Cualquier solición factible para el problema (2.80)-(2.82) es factible para el problema (2.76)-(2.79). En particular,  𝜆∗es solución factible de (2.76)-(2.79), y por tanto 𝑧𝑃 𝑀𝐿 =𝑐  𝜆∗=∑𝑚 𝑖=1 𝜋𝑖𝑏𝑖+∑𝐾 𝑘=1 𝜇𝑘≤𝑧𝑃 𝑀𝐿. Comprobación de optimalidad para el problema maestro lineal. Es necesario comprobar si (𝜋, 𝜇)es solución factible dual de (2.76)-(2.79). Esto implica comprobar que se cumple 𝑐𝑘𝑥−𝜋𝐴𝑘𝑥−𝜇𝑘≤0para cada 𝑥∈𝑋𝑘de cada 𝑘, de cada columna. En lugar de examinar cada punto de forma separada, se examinan todos los puntos 𝑋𝑘 implícitamente resolviendo el subproblema de optimización: 𝜁𝑘= m ax{(𝑐𝑘−𝜋𝐴𝑘)𝑥−𝜇𝑘∶𝑥∈𝑋𝑘}(2.84) Criterio de parada. Si 𝜁𝑘= 0 para 𝑘= 1,…, 𝐾, entonces la solución dual (𝜋, 𝜇)es la solución dual factible del problema maestro lineal y por tanto 𝑧𝑃 𝑀𝐿 ≤∑𝑚 𝑖=1 𝜋𝑖𝑏𝑖+ ∑𝐾 𝑘=1 𝜇𝑘. Como el valor de la solución factible primal  𝜆es igual al de su cota superior,  𝜆es óptimo para el problema maestro lineal. 40 trabajo fin de máster y proposiciones que aparecen en este capítulo se encuentran en el libro publicado por Ahuja,Magnanti yOrlin [3] en 1993. 3.2 Problemas de flujo con coste mínimo El modelo de flujo de mínimo coste es el más fundamental de todos los problemas de flujo en redes. El problema es fácil de enunciar: se desea determinar el envio de un servicio al menor coste a través de una red para satisfacer las demandas en ciertos nodos de los suministros disponibles en otros nodos. Sea 𝐺= (𝑁, 𝐴)una red dirigida definida por un conjunto 𝑁de nodos y un conjunto 𝐴de 𝑚arcos dirigidos. Cada arco (𝑖, 𝑗) ∈ 𝐴tiene asociado un coste 𝑐𝑖𝑗 que denota el coste por unidad de flujo en ese arco, luego el coste del flujo varía linealmente con la cantidad de flujo. A cada arco (𝑖, 𝑗) ∈ 𝐴también se asocia una capacidad 𝑢𝑖𝑗 que es la cantidad máxima de flujo en cada arco y un límite inferior 𝑙𝑖𝑗 que es la cantidad mínima de flujo que atraviesa ese arco. A cada nodo 𝑖∈𝑁se le asocia un número entero 𝑏(𝑖)que representa la oferta/demanda. Si 𝑏(𝑖)>0, el nodo 𝑖es un nodo de suministro; si 𝑏(𝑖)<0, el nodo 𝑖es un nodo demandante con una demanda −𝑏(𝑖); y si 𝑏(𝑖)=0, el nodo 𝑖es un nodo de transbordo. La variable de decisión en el problema de flujo de mínimo coste es el flujo 𝑥𝑖𝑗 en el arco (𝑖, 𝑗) ∈ 𝐴. El problema de flujo de mínimo coste es un modelo de optimización que se formula así: Minimizar ∑ (𝑖,𝑗)∈𝐴 𝑐𝑖𝑗𝑥𝑖𝑗 (3.1) s. a ∑ {𝑗∶(𝑖,𝑗)∈𝐴} 𝑥𝑖𝑗 −∑ {𝑗∶(𝑗,𝑖)∈𝐴} 𝑥𝑗𝑖 =𝑏(𝑖)para todo 𝑖∈𝑁, (3.2) 𝑙𝑖𝑗 ≤𝑥𝑖𝑗 ≤𝑢𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴, (3.3) donde ∑𝑛 𝑖=1 𝑏(𝑖)=0. El problema de flujo con mínimo coste en forma matricial es el siguiente: Minimizar 𝑐𝑥 (3.4) s. a 𝑥=𝑏, (3.5) 𝑙≤𝑥≤𝑢. (3.6) donde es una matriz 𝑛×𝑚llamada matriz de incidencia de los arcos en los nodos del problema de flujo con coste mínimo. La columna 𝑖𝑗 tiene un +1 en la 𝑖-ésima fila, un −1 en la 𝑗-ésima fila; el resto de sus filas son cero. 3. problema de flujo de mínimo coste: algoritmo simplex en redes 41 Las restricciones en (3.2) son restricciones de balance de flujo. El primer término de esta restricción representa el flujo de salida total del nodo (es decir, el flujo que emana del nodo) y el segundo término representa el flujo de entrada total del nodo (es decir, el flujo que ingresa al nodo). La restricción del balance de flujo establece que la salida menos la entrada debe ser igual a la oferta/demanda de ese nodo. Si el nodo es un nodo de suministro, su flujo de salida excede su flujo de entrada; si el nodo es un nodo de demanda, su entrada excede su salida; y si el nodo es un nodo de transbordo, su flujo de salida es igual a su flujo de entrada. El flujo también debe satisfacer el límite inferior y las restricciones de capacidad (3.3), que son las restricciones de límite de flujo. Los límites de flujo típicamente modelan capacidades físicas o restricciones impuestas en los rangos de operación de los flujos. En la mayoría de las aplicaciones, los límites inferiores de los flujos en arco son cero; por lo tanto, si no se establecen límites inferiores para ningún problema, se asumen que tienen valor cero. Se asume, en general, que los datos son números enteros, esto es el supuesto de integralidad. El cual, no es restrictivo para la mayoría de las aplicaciones porque siempre se puede transformar datos racionales en datos enteros multiplicándolos por una constante. Además, es necesario convertir números irracionales en números racionales para representarlos en un ordenador. Se enuncia a continuación varias variantes del problema (3.1) - (3.3) que aparecen posteriormente. Problema del camino más corto. Para este problema se desea encontrar el camino de mínimo coste (o longitud) que va desde un nodo fuente 𝑠a un nodo sumidero 𝑡, cada arco (𝑖, 𝑗) ∈ 𝐴tiene un coste (longitud) asociado 𝑐𝑖𝑗. Problema de circulación. El problema de circulación es un problema de flujo de coste mínimo con sólo nodos de transbordo; es decir, 𝑏(𝑖) = 0 para todo 𝑖∈𝑁. En este caso se desea encontrar un flujo factible que respete las cotas superior e inferior 𝑙𝑖𝑗 y𝑢𝑖𝑗 impuestas sobre las capacidades de los flujos de los arcos 𝑥𝑖𝑗. Dado que nunca se introduce ningún flujo exógeno en la red ni se extrae ningún flujo de ella, todo el flujo circula por la red. Se desea encontrar la circulación que tenga el coste mínimo. 42 trabajo fin de máster Problemas de flujo con costes convexos. En el problema de flujo con mínimo coste, el coste del flujo en cualquier arco varía linealmente con la cantidad de flujo. Los problemas de flujo de costes convexos tienen una estructura de costes más general: el coste es una función convexa de la cantidad de flujo. Problemas de flujo multiservicio. El problema de flujo con coste mínimo modela el flujo de un único servicio a través de una red. Los problemas de flujo multiservicio surgen cuando varios servicios utilizan la misma red subyacente. Los servicios pueden diferenciarse por sus características físicas o simplemente por sus pares de origen-destino. Los diferentes servicios tienen diferentes orígenes y destinos, y los servicios tienen restricciones de balance de flujo independientes en cada nodo. Sin embargo, comparten las capacidades de los arcos comunes que llevan varios servicios. De hecho, el problema esencial que aborda el problema del flujo multiservicio es la asignación de la capacidad de cada arco a los servicio individuales de una manera que minimice los costes generales del flujo. A continuación, se enumeran varios tipos de transformaciones en redes será útil para posteriores desarrollos. División de nodos. Esta transformación divide cada nodo 𝑖en dos nodos 𝑖′e𝑖′′ correspondientes a las funciones de entrada y salida del nodo. Esta transformación reemplaza cada arco original (𝑖, 𝑗)por un arco (𝑖′, 𝑗′′ )del mismo coste y capacidad. También agrega un arco (𝑖′′ , 𝑖′) de coste cero y capacidad infinita para cada 𝑖. La entrada del nodo 𝑖(es decir, el nodo 𝑖′′ ) recibe todo el flujo de entrada del nodo, la salida (es decir, el nodo 𝑖′) envía todo el flujo de salida del nodo, y el arco adicional (𝑖′′ , 𝑖′)lleva el flujo desde la entrada a la salida. La Figura 3.1 ilustra la red resultante cuando se lleva a cabo la transformación de división de nodos para todos los nodos de un red. Se define las ofertas/demandas de los nodos de la red transformada de acuerdo con los siguientes tres casos: 1. Si 𝑏(𝑖)>0, entonces 𝑏(𝑖′′ ) = 𝑏(𝑖)y𝑏(𝑖′)=0. 2. Si 𝑏(𝑖)<0, entonces 𝑏(𝑖′′ )=0y𝑏(𝑖′) = 𝑏(𝑖). 3. problema de flujo de mínimo coste: algoritmo simplex en redes 43 3. Si 𝑏(𝑖)=0, entonces 𝑏(𝑖′) = 𝑏(𝑖′′ )=0. (a) (b) Figura 3.1: División de nodos: red original(a); red transformada (b). Es fácil mostrar una correspondencia uno a uno entre un flujo en la red original y el flujo correspondiente en la red transformada; además, los flujos en ambas redes tienen el mismo coste. La transformación de división de nodos permite modelar numerosas aplicaciones en una variedad de dominios de problemas prácticos. Por ejemplo, se puede usar la transformación para manejar situaciones en las que los nodos y los arcos tienen capacidades y costes asociados. En estas situaciones, cada unidad de flujo que pasa por un nodo 𝑖incurre en un coste 𝑐𝑖y el flujo máximo que puede pasar a través del nodo es 𝑢𝑖. Puede reducirse este problema a la forma estándar de ”flujo de arco” del problema de flujo de red realizando la división de nodos y dejando que 𝑐𝑖y𝑢𝑖sean el coste y la capacidad respectivamente del arco (𝑖′′ , 𝑖′). 44 trabajo fin de máster Trabajar con costes reducidos Supóngase que cada nodo 𝑖∈𝑁tiene asociado un número 𝜋(𝑖)que es el potencial del nodo. Respecto a los potenciales del nodo 𝜋= (𝜋(1), 𝜋(2),…, 𝜋(𝑛)), se define el coste reducido 𝑐𝜋 𝑖𝑗 de un arco (𝑖, 𝑗)como: 𝑐𝜋 𝑖𝑗 =𝑐𝑖𝑗 −𝜋(𝑖) + 𝜋(𝑗)(2.3) Con frecuencia se trabaja con costes reducidos 𝑐𝜋 𝑖𝑗 en lugar de los costes reales 𝑐𝑖𝑗. Así, es importante entender la relación entre las funciones objetivo 𝑧(𝜋) = ∑(𝑖,𝑗)∈𝐴𝑐𝜋 𝑖𝑗𝑥𝑖𝑗 y𝑧(0) = ∑(𝑖,𝑗)∈𝐴𝑐𝑖𝑗𝑥𝑖𝑗. Supóngase que inicialmente 𝜋= 0, luego se incrementa el potencial de nodo del nodo 𝑘hasta 𝜋(𝑘). La definición (2.3) de costes reducidos supone que este cambio reduce el coste reducido de cada unidad de flujo que sale del nodo 𝑘 en 𝜋(𝑘)y aumenta el coste reducido de cada unidad de flujo que ingresa al nodo 𝑘en 𝜋(𝑘). Por tanto, la disminución total de la función objetivo es igual a 𝜋(𝑘)multiplicado por el flujo de salida del nodo 𝑘menos el flujo de entrada del nodo 𝑘. Por definición, la salida menos la entrada es igual a la oferta/demanda del nodo. En consecuencia, aumentar el potencial del nodo 𝑘en 𝜋(𝑘)disminuye el valor de la función objetivo en 𝜋(𝑘)𝑏(𝑘)unidades. Iterando para cada nodo se llega a: 𝑧(0) − 𝑧(𝜋) = ∑ 𝑖∈𝑁 𝜋(𝑖)𝑏(𝑖) = 𝜋𝑏. (3.7) Para un potencial de nodo dado 𝜋,𝜋𝑏 es una constante. Por tanto, un flujo que minimiza 𝑧(𝜋)también minimiza 𝜋(0). Propiedad 3.1.El problema de flujo de mínimo coste con costes en los arcos 𝑐𝑖𝑗 o𝑐𝜋 𝑖𝑗 tiene las mismas soluciones óptimas. Además, 𝑧(𝜋) = 𝑧(0) − 𝜋𝑏. A continuación, se estudia los efectos de trabajar con costes reducidos en los costes de ciclos y caminos. Sea 𝑊un ciclo dirigido en 𝐺. Entonces: ∑ (𝑖,𝑗)∈𝑊 𝑐𝜋 𝑖𝑗 =∑ (𝑖,𝑗)∈𝑊 (𝑐𝑖𝑗 −𝜋(𝑖) + 𝜋(𝑗)),(3.8) =∑ (𝑖,𝑗)∈𝑊 𝑐𝑖𝑗 +∑ (𝑖,𝑗)∈𝑊 (𝜋(𝑗) − 𝜋(𝑖)),(3.9) =∑ (𝑖,𝑗)∈𝑊 𝑐𝑖𝑗.(3.10) La última igualdad (3.10) resulta de que en cada ciclo dirijido 𝑊, la expresión ∑(𝑖,𝑗)∈𝑊(𝜋(𝑗)− 𝜋(𝑖)) suma cero para cada nodo 𝑖en el ciclo 𝑊porque 𝜋(𝑖)ocurre una vez con signo 3. problema de flujo de mínimo coste: algoritmo simplex en redes 45 positivo y otra con signo negativo. De forma similar, si 𝑃es un camino dirigido desde el nodo 𝑘hasta el nodo 𝑙, entonces: ∑ (𝑖,𝑗)∈𝑃 𝑐𝜋 𝑖𝑗 =∑ (𝑖,𝑗)∈𝑃 (𝑐𝑖𝑗 −𝜋(𝑖) + 𝜋(𝑗)),(3.11) =∑ (𝑖,𝑗)∈𝑃 𝑐𝑖𝑗 −∑ (𝑖,𝑗)∈𝑃 (𝜋(𝑖) + 𝜋(𝑗)),(3.12) =∑ (𝑖,𝑗)∈𝑃 𝑐𝑖𝑗 −𝜋(𝑘) + 𝜋(𝑙),(3.13) porque todos los 𝜋(⋅)correspondientes a los nodos en la ruta y distintos de los nodos terminales 𝑘e𝑙, se cancelan entre sí en la expresión ∑(𝑖,𝑗)∈𝑃(𝜋(𝑖) − 𝜋(𝑗)). Propiedad 3.2.1. Para todo ciclo dirigido 𝑊y para todo potencial 𝜋,∑(𝑖,𝑗)∈𝑊𝑐𝜋 𝑖𝑗 = ∑(𝑖,𝑗)∈𝑊𝑐𝑖𝑗. 2. Para cualquier camino dirigido 𝑃del nodo 𝑘hasta el nodo 𝑙y para cualquier potencial de nodo 𝜋,∑(𝑖,𝑗)∈𝑃𝑐𝜋 𝑖𝑗 =∑(𝑖,𝑗)∈𝑃𝑐𝑖𝑗 −𝜋(𝑘) + 𝜋(𝑙). Trabajando con redes residuales. Al diseñar, desarrollar e implementar algoritmos de flujo de red, a menudo es conveniente no medir el flujo en términos absolutos, sino en términos del incremento de flujo sobre alguna solución factible dada, típicamente, la solución en algún punto intermedio de un algoritmo. Hacerlo lleva a definir una nueva red auxiliar, conocida como red residual, que funciona como una ”red del flujo restante” para transportar el incremento de flujo. Se muestra que las formulaciones del problema en la red original y en la red residual son equivalentes, ambos dan una correspondencia biunívoca entre soluciones factibles a los dos problemas que preserva el valor del coste de las soluciones. El concepto de red residual se basa en la siguiente idea intuitiva. Supóngase que el arco (𝑖, 𝑗)lleva 𝑥0unidades de flujo. Entonces, se puede enviar 𝑢𝑖𝑗 −𝑥0 𝑖𝑗 unidades de flujo adicionales del nodo 𝑖al nodo 𝑗a lo largo del arco (𝑖, 𝑗). Obsérvese también que se puede enviar hasta 𝑥0unidades de flujo del nodo 𝑗al nodo 𝑖por el arco (𝑖, 𝑗), lo que equivale a cancelar el flujo existente en el arco. Mientras que enviar un flujo unitario del nodo 𝑖al nodo 𝑗en el arco (𝑖, 𝑗)aumenta el coste del flujo en 𝑐𝑖𝑗 unidades, enviar flujo desde el nodo 𝑗al nodo 𝑖en el mismo arco disminuye el coste del flujo en 𝑐𝑖𝑗 unidades (ya que se ahorra el coste en el que se solía incurrir al enviar el flujo del nodo 𝑖al nodo 𝑗). 46 trabajo fin de máster Usando estas ideas, se define la red residual con respecto a un flujo dado 𝑥0de la siguiente manera. Se reemplaza cada arco (𝑖, 𝑗)en la red original por dos arcos, (𝑖, 𝑗) y(𝑗, 𝑖): el arco (𝑖, 𝑗)tiene un coste 𝑐𝑖𝑗 y una capacidad residual 𝑟𝑖𝑗 =𝑢𝑖𝑗 −𝑥0 𝑖𝑗, y el arco (𝑗, 𝑖)tiene coste −𝑐𝑖𝑗 y capacidad residual 𝑟𝑗𝑖 =𝑥0 𝑖𝑗 (ver Figura 3.2). La red residual consiste sólo en los arcos con una capacidad residual positiva. Se usa la notación 𝐺(𝑥0) para representar la red residual correspondiente al flujo 𝑥0. Figura 3.2: Construcción de la red residual 𝐺(𝑥0). Ahora, se muestra que todo flujo 𝑥en la red 𝐺corresponde a un flujo 𝑥′en la red residual 𝐺(𝑥0). Se define el flujo 𝑥′≥0de la siguiente manera: 𝑥′ 𝑖𝑗 −𝑥′ 𝑗𝑖 =𝑥𝑖𝑗 −𝑥0 𝑖𝑗,(3.14) y 𝑥′ 𝑖𝑗𝑥′ 𝑗𝑖 = 0.(3.15) La condición (3.15) supone que 𝑥′ 𝑖𝑗 y𝑥′ 𝑗𝑖 no pueden ser ambos positivos al mismo tiempo. Si 𝑥𝑖𝑗 ≥𝑥0 𝑖𝑗, se establece 𝑥′ 𝑖𝑗 = (𝑥𝑖𝑗 −𝑥0 𝑖𝑗)y𝑥′ 𝑗𝑖 = 0. Obsérvese que si 𝑥𝑖𝑗 < 𝑢𝑖𝑗, entonces 𝑥′ 𝑖𝑗 ≤𝑢𝑖𝑗 −𝑥0 𝑖𝑗 =𝑟𝑖𝑗. Por lo tanto, el flujo 𝑥′ 𝑖𝑗 satisface las restricciones de límite de flujo. De manera similar, si 𝑥𝑖𝑗 < 𝑥0 𝑖𝑗, se establece 𝑥′ 𝑖𝑗 = 0 y𝑥′ 𝑗𝑖 =𝑥0 𝑖𝑗 −𝑥𝑖𝑗. Si 0≤𝑥′ 𝑗𝑖 ≤𝑥0 𝑖𝑗 =𝑟𝑗𝑖, entonces el flujo 𝑥′ 𝑗𝑖 también satisface las restricciones de límite de flujo. Estas observaciones muestran que si 𝑥es un flujo factible en 𝐺, su correspondiente flujo 𝑥′es un flujo factible en 𝐺(𝑥0). A continuación, se establece una relación entre el coste de un flujo 𝑥en 𝐺y el coste del flujo correspondiente 𝑥′en 𝐺(𝑥0). Sea 𝑐′los costes del arco en la red residual. Luego, para cada arco (𝑖, 𝑗) ∈ 𝐴,𝑐′ 𝑖𝑗 =𝑐𝑖𝑗 and 𝑐′ 𝑗𝑖 = −𝑐𝑖𝑗. Para un flujo 𝑥𝑖𝑗 en el arco (𝑖, 𝑗)en la red original 𝐺, el coste del flujo en el par de arcos (𝑖, 𝑗)y(𝑗, 𝑖)en la red residual 𝐺(𝑥0)es 𝑐′ 𝑖𝑗𝑥′ 𝑖𝑗 +𝑐′ 𝑗𝑖𝑥′ 𝑗𝑖 =𝑐′ 𝑖𝑗(𝑥′ 𝑖𝑗 −𝑥′ 𝑗𝑖) = 𝑐𝑖𝑗𝑥𝑖𝑗 −𝑐𝑖𝑗𝑥0 𝑖𝑗; la última igualdad se sigue de (3.14). Así se tiene que: 𝑐′𝑥′=𝑐𝑥 −𝑐𝑥0.(3.16) 3. problema de flujo de mínimo coste: algoritmo simplex en redes 47 De manera similar, se puede mostrar el resultado inverso de que si 𝑥′es un flujo factible en la red residual 𝐺(𝑥0), la solución dada por 𝑥𝑖𝑗 = (𝑥′ 𝑖𝑗 −𝑥′ 𝑗𝑖) + 𝑥0 𝑖𝑗 es un flujo factible en 𝐺. Además, los costes de estos dos flujos está relacionado por la igualdad 𝑐𝑥 =𝑐′𝑥′+𝑐𝑥0. Propiedad 3.3.Un flujo 𝑥es un flujo factible en la red 𝐺si y sólo si su correspondiente flujo 𝑥′, definido por 𝑥′ 𝑖𝑗 −𝑥′ 𝑗𝑖 =𝑥𝑖𝑗 +𝑥0 𝑖𝑗, es factible en la red residual 𝐺(𝑥0). Además, 𝑐𝑥 =𝑐′𝑥′+𝑐𝑥0. Una consecuencia importante de la Propiedad 3.3 es la flexibilidad que proporciona. En lugar de trabajar con la red 𝐺original, permite trabajar con la red residual 𝐺(𝑥0) para algún 𝑥0: una vez que se ha determinado una solución óptima en la red residual, se la puede convertir inmediatamente en una solución óptima en la red original. 3.3 Soluciones de árbol de expansión y soluciones sin ciclos Para cualquier solución factible, 𝑥, un arco (𝑖, 𝑗)es libre si 0< 𝑥𝑖𝑗 < 𝑢𝑖𝑗 y es un arco restringido si 𝑥𝑖𝑗 = 0 o𝑥𝑖𝑗 =𝑢𝑖𝑗. Puede aumentarse y disminuirse el flujo en un arco libre respetando las restricciones de flujo del arco. En la cota inferior de flujo de un arco restringido (𝑖, 𝑗)(𝑥𝑖𝑗 = 0) sólo se puede aumentar el flujo. De manera similar, en la cota superior de flujo de un arco restringido (𝑖, 𝑗)(𝑥𝑖𝑗 =𝑢𝑖𝑗) sólo se puede disminuir el flujo. La solución 𝑥es una solución sin ciclos si la red no contiene ningún ciclo formado únicamente por arcos libres. En una solución sin ciclos, al aumentar el flujo en un ciclo de aumento (ver Sección 5.3.1) sólo se puede hacer en una dirección, puesto que en el ciclo siempre habrá algún arco que impedirá aumentar o disminuir el flujo de ese arco. Una solución 𝑥es factible y es un árbol de expansión asociado a la red (solución en árbol de expansión) si cada arco que no es del árbol es un arco restringido. En una solución de árbol de expansión, los arcos que no son del árbol pueden ser libres o restringidos. Con frecuencia, cuando se hace referencia a una solución de árbol de expansión, no se identifica explícitamente el árbol asociado; más bien, se entenderá por el contexto. Los problemas de flujo de mínimo coste siempre tienen soluciones óptimas de árbol de expansión y sin ciclo. El algoritmo simplex en redes explota este resultado al restringir su búsqueda a sólo soluciones de árbol de expansión. 48 trabajo fin de máster (a) (b) Figura 3.3: a) Solución factible; b) solución tras aumentar una cantidad 𝜃de flujo en el ciclo. Supóngase que los arcos no están capacitados (𝑢𝑖𝑗 = ∞ para cada (𝑖, 𝑗) ∈ 𝐴). La red que se muestra en la Figura 3.3 contiene flujo positivo alrededor de un ciclo. Se define la orientación del ciclo como la misma que la del arco (4,5). Se aumenta 𝜃unidades de flujo a lo largo del ciclo en la dirección de su orientación. Como se muestra en la Figura 3.3, esto aumenta el flujo en los arcos a lo largo de la orientación del ciclo (arcos hacia adelante) en 𝜃unidades y disminuye el flujo en los arcos opuestos a la orientación del ciclo (arcos hacia atrás) en 𝜃unidades. El coste incremental por unidad para este cambio de flujo es la suma de los costes de los arcos hacia adelante menos la suma de los costes de los arcos hacia atrás en el ciclo, es decir, Aumento de coste por unidad de flujo: Δ = 2 + 1 + 3 − 4 − 3 = −1.(3.17) Dado que aumentar el flujo en el ciclo disminuye el coste, se establece 𝜃lo más grande posible mientras se preserva la no negatividad de todos los flujos en los arcos. Se debe satisfacer las desigualdades 3 − 𝜃≥0y4 − 𝜃≥0, y por lo tanto se establece 𝜃= 3. En la nueva solución (en 𝜃= 3), algún arco en el ciclo tiene un flujo con valor cero, y además, esta solución es estrictamente menor que el valor de la solución inicial. Si se cambia 𝑐12 de 2a5, el coste unitario del ciclo es Δ=2. En consecuencia, para mejorar al máximo el coste, se disminuye 𝜃tanto como sea posible (satisfaciendo las restricciones 5 + 𝜃≥0,2 + 𝜃≥0,y4 + 𝜃≥0,o𝜃≥−2) y se encuentra nuevamente la solución de menor coste con al menos un arco con flujo nulo. Puesto de otra manera: 3. problema de flujo de mínimo coste: algoritmo simplex en redes 49 para preservar todos los flujos en arco como no negativos, debe seleccionarse 𝜃en el intervalo −2 ≤𝜃≤3. Dado que la función objetivo depende linealmente de 𝜃, se optimiza seleccionando 𝜃= 3 o𝜃= −2, en este punto un arco en el ciclo tiene valor nulo. Se puede extender esto de la siguientes formas: 1. Si el incremento de coste Δ=0, no se obtiene mejora con las soluciones en el intervalo −2 ≤𝜃≤3y, por lo tanto, vuelve a elegirse una solución tan buena como la original, pero con el flujo de al menos un arco en el ciclo con valor cero. 2. Si se limita superiormente el flujo (por ejemplo no más de 6 unidades en cada arco), el rango de flujos que conserva la viabilidad (las restricciones de balance de flujo, las cotas inferior y superior de los flujos) es nuevamente un intervalo, en este caso: −2 ≤𝜃≤1, y puede encontrarse una solución tan buena como la original eligiendo 𝜃= −2 o𝜃= 1. Con estos valores de 𝜃, la solución no tiene ciclos; es decir, algún arco en el ciclo tiene un flujo con valor cero. Teorema 3.1 (Propiedad sin ciclos). Si la función objetivo de un problema de flujo de mínimo coste está acotada inferiormente dentro de la región factible, entonces el problema siempre tiene una solución óptima sin ciclos. Los arcos libres en una solución sin ciclos definen un bosque (una colección de árboles cuyos nodos están separados). Si este bosque es un árbol de expansión, la solución sin ciclo ya es una solución de árbol de expansión. La Figura 3.4 (c) corresponde a un árbol de expansión de una solución sin ciclos. Podría ser posible (y a menudo lo es) formar de varias formas un árbol de expansión a partir del conjunto de arcos. Agregar el arco (3,4) en lugar del arco (2,4) o(3,5) produciría otra solución de árbol de expansión. Por lo tanto, una solución sin ciclo dada puede corresponder a varios árboles de expansión. Sin embargo, dado que se asume que la red subyacente es conexa, siempre pueden agregarse algunos arcos restringidos a los arcos libres de una solución sin ciclos para producir un árbol de expansión, por lo que se establece el siguiente resultado: Teorema 3.2 (propiedad del árbol de expansión). Si la función objetivo de un problema de flujo de mínimo coste está acotado inferiormente dentro de la región factible, el problema siempre tiene una solución óptima de árbol de expansión. Una solución de árbol de expansión divide el conjunto de arcos 𝐴en tres subconjuntos: (1) 𝑇, los arcos en el árbol de expansión; (2) 𝐿, los arcos que no son del árbol cuyo flujo está restringido a cero; y (3) 𝑈, los arcos que no son del árbol cuyo flujo está 56 trabajo fin de máster caso en el que el nodo 𝑖es el nodo hoja es similar), el arco (𝑖, 𝑗)debe transportar −𝑏(𝑗) unidades de flujo. El Algoritmo 2 ofrece una descripción formal de este procedimiento. 𝑏′(𝑖)←𝑏(𝑖),para todo 𝑖∈𝑁; 𝑇′ ←𝑇; para (𝑖, 𝑗) ∈ 𝐿hacer 𝑥𝑖𝑗 ←0 fin 𝑗←𝑟𝑎𝑚𝑎(𝑖); mientras 𝑇′≠{1} hacer seleccionar un nodo hoja 𝑗(que no sea el nodo 1) en el subárbol 𝑇′; 𝑖←𝑝𝑟𝑒𝑑(𝑗); si (𝑖, 𝑗) ∈ 𝑇′entonces 𝑥𝑖𝑗 ←−𝑏′(𝑗) en otro caso 𝑥𝑖𝑗 ←𝑏′(𝑗) fin 𝑏′(𝑖)←𝑏′(𝑖) + 𝑏′(𝑗); eliminar el nodo 𝑗y el arco indicente al este nodo y parte de 𝑇′; fin Algoritmo 2: Procedimiento de cálculo de flujos Este método para calcular los flujos asume que el problema de flujo de mínimo coste no está capacitado. Para la versión capacitada del problema, se agrega el siguiente enunciado inmediatamente después de la primera condición (𝑏′(𝑖) ∶= 𝑏(𝑖)para todo 𝑖∈𝑁) en el procedimiento de cálculo de flujos. para (𝑖, 𝑗) ∈ 𝑈hacer 𝑥𝑖𝑗 ←𝑢𝑖𝑗; 𝑏′(𝑖)←𝑏′(𝑖) − 𝑢𝑖𝑗; 𝑏′(𝑗)←𝑏′(𝑗) + 𝑢𝑖𝑗; fin La inicialización de los flujos y la modificación de las ofertas/demandas 𝑏(𝑖)y𝑏(𝑗)para los arcos (𝑖, 𝑗)en 𝑈requiere un tiempo (𝑚). Ignorando el tiempo para seleccionar los nodos hoja de 𝑇cada iteración requiere tiempo (1), lo que da como resultado un total de tiempo (𝑛). Una forma de identificar los nodos hoja en 𝑇es seleccionar nodos en el orden inverso al de los índices rama. Cada nodo aparece antes que sus descendientes en el recorrido de la rama. Este recorrido inverso examina cada nodo sólo después de 3. problema de flujo de mínimo coste: algoritmo simplex en redes 57 haber examinado todos los descendientes del nodo. Así, para el problema de flujo de mínimo coste sin capacidades, el procedimiento de cálculo de flujo se ejecuta en un tiempo (𝑚). Para la versión con capacidades del problema, el procedimiento también requiere un tiempo (𝑚). El problema de flujo de mínimo coste siempre tiene una solución óptima de árbol de expansión. El flujo asociado con este árbol es un flujo óptimo y puede determinarse usando el procedimiento de cálculo de flujo. Este procedimiento implica que si las capacidades de todos los arcos y las ofertas/demandas de todos los nodos son enteras, los flujos de arco también son enteros (el procedimiento solo realiza sumas y restas). Teorema 3.5 (Propiedad de integralidad primaria). Si las capacidades de todos los arcos y las ofertas/demandas de todos los nodos son enteras, el problema de flujo de mínimo coste tiene siempre un flujo óptimo entero. Cada estructura de árbol de expansión (𝑇 , 𝐿, 𝑈)define un flujo 𝑥único. Si este flujo satisface las restricciones de flujo 0≤𝑥𝑖𝑗 ≤𝑢𝑖𝑗 para cada arco (𝑖, 𝑗) ∈ 𝐴, la estructura del árbol de expansión es factible; de lo contrario, no es factible. El árbol de expansión 𝑇es degenerado si 𝑥𝑖𝑗 = 0 o𝑥𝑖𝑗 =𝑢𝑖𝑗 para algún arco (𝑖, 𝑗) ∈ 𝑇, y no degenerado en caso contrario. Para cada arco (𝑖, 𝑗)en un árbol de expansión no degenerado se cumple 0< 𝑥𝑖𝑗 < 𝑢𝑖𝑗. 3.6 Algoritmo simplex en redes El algoritmo simplex en redes mantiene una estructura de árbol de expansión factible en cada iteración y transforma esta estructura sucesivamente en una estructura de árbol de expansión mejorada hasta que se vuelve óptima. El Algoritmo 4 especifica los pasos esenciales del método. 58 trabajo fin de máster determinar una estructura de árbol factible inicial (𝑇 , 𝐿, 𝑈); sea 𝑥el flujo y 𝜋el potencial de nodo asociado a la estructura del árbol; mientras algunos arcos no pertenecientes al árbol violen la condición de optimalidad hacer seleccionar un arco de entrada (𝑘, 𝑙)que viola la condición de optimalidad; añade el arco (𝑘, 𝑙)al árbol y determinar el arco saliente (𝑝, 𝑞); realizar una actualización de árbol y actualizar las soluciones 𝑥y𝜋; fin Algoritmo 3: Algoritmo simplex en redes 3.6.1 Obtención de una estructura de árbol de expansión inicial Se asume que para cada nodo 𝑗∈𝑁− {1}, la red contiene los arcos (1, 𝑗)y(𝑗, 1), con costes y capacidades suficientemente grandes. El árbol inicial 𝑇se construye de la siguiente manera. Se examina cada nodo 𝑗, distinto del nodo 1, uno por uno. Si 𝑏(𝑗)≥0, se incluye el arco (1, 𝑗)en 𝑇con un valor de flujo 𝑏(𝑗). Si 𝑏(𝑗)<0, se incluye el arco (𝑗, 1) en 𝑇con un valor de flujo de −𝑏(𝑗). El conjunto 𝐿consta de los arcos restantes y el conjunto 𝑈está vacío. Puede calcularse fácilmente los potenciales de nodo para este árbol usando las ecuaciones 𝑐𝑖𝑗 −𝜋(𝑖)+𝜋(𝑗)=0para todo (𝑖, 𝑗) ∈ 𝑇. Téngase en cuenta que 𝜋(1) = 0. Si la red no contiene los arcos (1, 𝑗)y(𝑗, 1) para cada nodo 𝑗∈𝑁− {1} (o no se desea agregar estos arcos por alguna razón), puede construirse una estructura de árbol de expansión inicial estableciendo primero un flujo factible en la red resolviendo un problema de flujo máximo y luego convirtiendo esta solución en una solución de árbol de expansión. 3.6.2 Condiciones de optimalidad Sea (𝑇 , 𝐿, 𝑈)una estructura de árbol de expansión factible del problema de flujo de mínimo coste, y sea 𝜋los correspondientes potenciales de nodo. Para determinar si la estructura de árbol de expansión es óptima, se verifica si cumple las siguientes 3. problema de flujo de mínimo coste: algoritmo simplex en redes 59 condiciones: 𝑐𝜋 𝑖𝑗 ≥0para todo arco (𝑖, 𝑗) ∈ 𝐿, (3.24) 𝑐𝜋 𝑖𝑗 ≤0para todo arco (𝑖, 𝑗) ∈ 𝑈. (3.25) Si la estructura del árbol de expansión satisface estas condiciones, es óptima y el algoritmo termina. De lo contrario, el algoritmo selecciona un arco no perteneciente al árbol que viola la condición de optimalidad para ser introducido en el árbol. Dos tipos de arcos son seleccionables para ingresar al árbol: 1. Cualquier arco (𝑖, 𝑗) ∈ 𝐿con 𝑐𝜋 𝑖𝑗 <0, 2. Cualquier arco (𝑖, 𝑗) ∈ 𝑈con 𝑐𝜋 𝑖𝑗 >0. Para cualquier arco seleccionable (𝑖, 𝑗),|𝑐𝜋 𝑖𝑗|es su violación. El algoritmo simplex en redes puede seleccionar cualquier arco elegible para ingresar en el árbol y aún así terminaría de manera finita (con algunas disposiciones para lidiar con la degeneración). Sin embargo, diferentes reglas para seleccionar el arco de entrada producen algoritmos con diferente comportamiento empírico y teórico. Son posibles muchas reglas diferentes, llamadas reglas de pivote, para elegir el arco de entrada. Las siguientes reglas son las más ampliamente adoptadas. Regla de pivote de Dantzig. Esta regla selecciona en cada iteracción el arco con mayor violación puesto que causa la máxima disminución en la función objetivo por unidad de flujo. Si el aumento promedio en el flujo del arco seleccionado fuera el mismo para todos los arcos entonces la introducción de este arco en el árbol de expansión causa la máxima disminución por pivote. Esta elección del arco de entrada tiende a producir disminuciones relativamente grandes en la función objetivo por iteración y, como resultado, el algoritmo realiza menos iteraciones que otras opciones para la regla de pivote. Sin embargo, esta regla tiene un gran inconveniente: el algoritmo debe considerar cada arco no perteneciente al árbol para identificar el arco con la máxima violación y hacerlo requiere mucho tiempo. Por lo tanto, aunque este algoritmo generalmente realiza menos iteraciones que otras implementaciones, el tiempo de ejecución del algoritmo no es atractivo. Primera regla de pivote de arco elegible. Para implementar esta regla, se escanea la lista de arcos secuencialmente y se selecciona el primer arco apropiado para ingresar en el árbol. En una versión popular de esta regla, se examina la lista de la forma en la que a continuación se ejemplifica: en una iteración, si se encuentra que el quinto arco en la lista de arcos es el primer arco apropiado, en la siguiente iteración 60 trabajo fin de máster comienza a escanearse la lista de arcos desde el sexto arco. Si se llega al final de la lista de arcos mientras se está realizando alguna iteración, se continua examinando la lista de arcos desde el principio. Una característica interesante de esta regla de pivote es que identifica rápidamente el arco de entrada. Esta regla de pivote tiene un inconveniente de contrapeso: con ella, el algoritmo generalmente realiza más iteraciones que con otras reglas de pivote porque cada operación de pivote produce una disminución relativamente pequeña en el valor de la función objetivo. El tiempo de ejecución del algoritmo no es muy atractivo, aunque la regla produce una implementación más eficiente que la regla de pivote de Dantzig. Lista de candidatos según la regla del pivote. El algoritmo selecciona el arco de entrada usando un procedimiento de dos fases que consta de iteraciones mayores e iteraciones menores. En una iteración mayor se construye una lista de candidatos de arcos elegibles. Una vez construida esta lista, se realiza una serie de iteraciones menores; en cada una de estas iteraciones, se selecciona un arco elegible de la lista de candidatos que causa la máxima disminución por unidad de flujo en la función objetivo. En una iteración mayor, se construye la lista de candidatos de la siguiente manera. Primero se examinan los arcos que parten del nodo 1 y se agregan arcos elegibles a la lista de candidatos. Este proceso se repite para los nodos 2,3,…,hasta que la lista haya alcanzado el máximo tamaño permitido o se hayan examinado todos los nodos. La siguiente iteración principal comienza con el nodo donde terminó la iteración principal anterior. El método de la lista de candidatos ofrece una flexibilidad considerable para ajustarse con precisión a la clase de problema. Al establecer el tamaño máximo permitido de la lista de candidatos de manera adecuada y al especificar el número de iteraciones menores que se realizarán dentro de una iteración principal, puede obtenerse numerosas reglas de pivote diferentes. A continuación, se estudia el proceso de elección del arco que abandona la estructura del árbol de expansión en cada paso del algoritmo simplex en redes. 3.6.3 Arcos que abandonan el árbol Se selecciona el arco de entrada (𝑘, 𝑙), la inclusión de este arco al árbol 𝑇crea un ciclo pivote 𝑊. El ciclo pivote consiste en la única ruta en el árbol 𝑇que va desde 3. problema de flujo de mínimo coste: algoritmo simplex en redes 61 el nodo 𝑘al nodo 𝑙, junto con el arco (𝑘, 𝑙). Se define la orientación del ciclo como la misma que (𝑘, 𝑙)si (𝑘, 𝑙) ∈ 𝐿y opuesta a la de (𝑘, 𝑙)si (𝑘, 𝑙) ∈ 𝑈. Sean 𝑊y𝑊 los conjuntos de arcos hacia adelante (arcos de 𝑊con su misma orientación) y arcos hacia atrás (arcos de 𝑊con la orientación contraria a 𝑊) en el ciclo pivote. El envío de flujo adicional por el ciclo pivote 𝑊en la dirección de su orientación disminuye estrictamente el coste de la solución actual en una tasa de |𝑐𝜋 𝑘𝑙|por unidad. Se aumenta el flujo tanto como sea posible hasta que uno de los arcos en el ciclo pivote alcance su cota inferior o superior. Aumentar el flujo a lo largo de 𝑊aumenta el flujo en los arcos hacia adelante y disminuye el flujo en los arcos hacia atrás. En consecuencia, el cambio de flujo máximo 𝛿𝑖𝑗 en un arco (𝑖, 𝑗) ∈ 𝑊que satisface las restricciones de límite de flujo es: 𝛿𝑖𝑗 ={𝑢𝑖𝑗 −𝑥𝑖𝑗 si (𝑖, 𝑗) ∈ 𝑊 , 𝑥𝑖𝑗 si (𝑖, 𝑗) ∈ 𝑊 . (3.26) Puede aumentarse 𝛿= m ın{𝛿𝑖𝑗 ∶ (𝑖, 𝑗) ∈ 𝑊}unidades de flujo a lo largo de 𝑊para mantener la viabilidad. Cualquier arco (𝑖, 𝑗) ∈ 𝑊que defina 𝛿(𝛿=𝛿𝑖𝑗)es un arco de bloqueo. Se aumenta 𝛿unidades de flujo y se selecciona un arco (𝑝, 𝑞)con 𝛿𝑝𝑞 =𝛿 como arco saliente. Una iteración pivote es una iteración no degenerada si 𝛿 > 0y es una iteración degenerada si 𝛿= 0. Una iteración degenerada ocurre sólo si 𝑇es un árbol de expansión degenerado. Si dos arcos tienen el mismo valor de 𝛿, el siguiente árbol de expansión se degenerará. El paso crucial para identificar el arco saliente es identificar el ciclo de pivote. Si 𝑃(𝑖) denota la única ruta en el árbol desde un nodo 𝑖al nodo raíz, este ciclo consta de los arcos {(𝑘, 𝑙)} ∪ 𝑃(𝑘) ∪ 𝑃(𝑙) − (𝑃(𝑘) ∩ 𝑃(𝑙)). En otras palabras, 𝑊consiste en el arco (𝑘, 𝑙)y las porciones disjuntas de 𝑃(𝑘)y𝑃(𝑙). Los índices predecesor permiten identificar el ciclo que sigue. Primero, se designan todos los nodos de la red como ”sin marcar”. Luego se comienza en el nodo 𝑘y, usando los índices predecesor, se trazan la ruta desde este nodo hasta la raíz y se marcan todos los nodos en esta ruta. Luego partiendo del nodo 𝑙se trazan los índices predecesor hasta que se encuentra un nodo marcado 𝑤. El nodo 𝑤es el primer ancestro común de los nodos 𝑘e𝑙; se trata del vértice esquina del ciclo 𝑊. El ciclo 𝑊contiene las porciones de las rutas 𝑃(𝑘)y𝑃(𝑙)hasta el nodo 𝑤, junto con el arco (𝑘, 𝑙). Este método identifica el tiempo del ciclo 𝑊en (𝑛), y por lo que es eficiente. Sin embargo, tiene el inconveniente de retroceder a lo largo de esos arcos de 𝑃(𝑘)que no están en 𝑊. Si el ciclo de pivote está lejos de su raíz, entonces rastrear los nodos hasta la raíz será ineficaz. Idealmente, se desea identificar el ciclo 𝑊en el tiempo proporcional a |𝑊|. El uso simultáneo de los índices profundidad y predecesor permite alcanzar este objetivo (veáse Algoritmo 4). 62 trabajo fin de máster 𝑖←𝑘y𝑗←𝑙; mientras 𝑖≠𝑗hacer si 𝑝𝑟𝑜𝑓(𝑖)> 𝑝𝑟𝑜𝑓(𝑗)entonces 𝑖←𝑝𝑟𝑒𝑑(𝑖) en otro caso si 𝑝𝑟𝑜𝑓(𝑗)> 𝑝𝑟𝑜𝑓(𝑖)entonces 𝑗←𝑝𝑟𝑒𝑑(𝑗) en otro caso 𝑖←𝑝𝑟𝑒𝑑(𝑖)y𝑗←𝑝𝑟𝑒𝑑(𝑗) fin fin fin Algoritmo 4: Procedimiento para identificar el ciclo pivote. Este método escanea los arcos en el ciclo de pivote 𝑊dos veces. Durante la primera exploración, se identifica el vértice esquina del ciclo y también se identifica el máximo flujo posible que puede aumentarse a lo largo de 𝑊. En la segunda exploración, se aumenta el flujo. Toda la operación de cambio de flujo requiere en el peor de los casos un tiempo (𝑛), pero normalmente examina solo un pequeño subconjunto de nodos (y arcos). 3.6.4 Actualización del árbol Cuando para un arco entrante dado (𝑘, 𝑙)el algoritmo simplex en redes ha determinado un arco saliente (𝑝, 𝑞)actualiza la estructura del árbol. Si el arco de salida es el mismo que el arco de entrada, lo que sucede cuando 𝛿=𝛿𝑘𝑙 =𝑢𝑘𝑙, el árbol no cambia. En este caso, el arco (𝑘, 𝑙)se mueve del conjunto 𝐿al conjunto 𝑈, o viceversa. Si el arco saliente difiere del arco entrante, el algoritmo debe realizar cambios más extensos. En este caso, el arco (𝑝, 𝑞)se convierte en un arco no perteneciente al árbol dependiendo de si 𝑥𝑝𝑞 = 0 o𝑥𝑝𝑞 =𝑢𝑝𝑞 (flujo actualizado) en su cota inferior o superior. Agregar el arco (𝑘, 𝑙)al árbol de expansión actual y eliminar el arco (𝑝, 𝑞)crea un nuevo árbol de expansión. Para el nuevo árbol de expansión, los potenciales de los nodos también cambian; pueden actualizarse de la siguiente manera. La eliminación del arco (𝑝, 𝑞)del árbol anterior divide el conjunto de nodos en dos subárboles. Uno, 𝑇1, que contiene el nodo raíz, y el otro, 𝑇2, que no contiene el nodo raíz. Téngase en cuenta que el subárbol 𝑇2 parte del nodo 𝑝o del nodo 𝑞. El arco (𝑘, 𝑙)tiene un extremo en 𝑇1y el otro en 𝑇2. 3. problema de flujo de mínimo coste: algoritmo simplex en redes 63 Como es fácil de verificar, las condiciones 𝜋(1) = 0 y𝑐𝑖𝑗 −𝜋(𝑖) + 𝜋(𝑗)=0implican para todos los arcos en el nuevo árbol que los potenciales de los nodos en el subárbol 𝑇1permanecen sin cambios, y los potenciales de los nodos en el subárbol 𝑇2cambian en una cantidad constante. Si 𝑘∈𝑇1e𝑙∈𝑇2, todos los potenciales de nodo en 𝑇2au- mentan en −𝑐𝜋 𝑘𝑙; si 𝑙∈𝑇1y𝑘∈𝑇2, aumentan en la cantidad 𝑐𝜋 𝑘𝑙. Usando los índices de profundidad y de rama, el método descrito en el Algoritmo 5 actualiza los potenciales de nodo rápidamente. si 𝑞∈𝑇2entonces 𝑦←𝑞 en otro caso 𝑦←𝑝 fin si 𝑘∈𝑇1entonces 𝑖𝑛𝑡𝑒𝑟𝑐𝑎𝑚𝑏𝑖𝑜 ←−𝑐𝜋 𝑘𝑙 en otro caso 𝑖𝑛𝑡𝑒𝑟𝑐𝑎𝑚𝑏𝑖𝑜 ←𝑐𝜋 𝑘𝑙 fin 𝜋(𝑦)←𝜋(𝑦) + 𝑖𝑛𝑡𝑒𝑟𝑐𝑎𝑚𝑏𝑖𝑜; 𝑧←𝑟𝑎𝑚𝑎(𝑦); mientras 𝑝𝑟𝑜𝑓(𝑧)> 𝑝𝑟𝑜𝑓(𝑦)hacer 𝜋(𝑧) = 𝜋(𝑧) + 𝑖𝑛𝑡𝑒𝑟𝑐𝑎𝑚𝑚𝑏𝑖𝑜; 𝑧←𝑟𝑎𝑚𝑎(𝑧); fin Algoritmo 5: Actualización del potencial de nodo en una operación pivote. El último paso en la actualización del árbol es volver a calcular los distintos índices del árbol. Este paso es bastante complicado pero es posible actualizar los índices del árbol en tiempo (𝑛). De hecho, el tiempo necesario para actualizar los índices del árbol es (|𝑊|+ m ın{|𝑇1|,|𝑇2|}), que suele ser mucho menor que 𝑛. 3.6.5 Terminación El algoritmo simplex en redes se mueve de una estructura de árbol de expansión factible a otra hasta que obtiene una estructura de árbol de expansión que satisface la condición de optimalidad (3.18)-(3.20). Si cada operación de pivote en el algoritmo no es degenerada, se puede mostrar que el algoritmo es finito. |𝑐𝜋 𝑘𝑙|es la disminución neta en el coste por flujo unitario enviado alrededor del ciclo de pivote 𝑊. En un pivote no degenerado (para el cual 𝛿 > 0), el coste de la nueva estructura de árbol de expansión 64 trabajo fin de máster es 𝛿|𝑐𝜋 𝑘𝑙|unidades menos que el coste de la estructura de árbol de expansión anterior. Dado que cualquier red tiene un número finito de estructuras de árbol de expansión y cada estructura de árbol de expansión tiene un coste asociado único, el algoritmo simplex en redes encuentra cualquier estructura de árbol de expansión como máximo una vez y, por lo tanto, terminará de manera finita. Sin embargo, los pivotes degenerados plantean una dificultad teórica: es posible que el algoritmo no termine de forma finita a menos que se realicen pivotes con cuidado. En la siguiente sección se discute una implementación especial, llamada implementación de árbol de expansión altamente factible, que garantiza la convergencia finita del algoritmo simplex en redes incluso para problemas que están degenerados. Se hace uso del ejemplo de la Figura 3.6 (a) para ilustrar el algoritmo simplex en redes. La Figura 3.6 (b) muestra una solución factible de árbol de expansión para el problema. Para esta solución, 𝑇= {(1,2),(1,3),(2,4),(2,5),(5,6)},𝐿= {(2,3),(5,4)} y𝑈= {(3,5),(4,6)}. En esta solución, el arco (3,5) tiene una violación positiva, que es una unidad. Se introduce este arco en el árbol creando un ciclo cuyo vértice es el nodo 1. Dado que el arco (3,5) está en su cota superior, la orientación del ciclo es opuesta a la del arco (3,5). Los arcos (1,2) y(2,5) son arcos hacia adelante en el ciclo y los arcos (3,5) y(1,3) son arcos hacia atrás. El aumento máximo de flujo permitido por los arcos (3,5),(1,3),(1,2) y(2,5) es, respectivamente, 3, 3, 2 y 1 unidades. En consecuencia, 𝛿= 1 y se aumenta una unidad de flujo a lo largo del ciclo. El aumento incrementa el flujo en los arcos (1,2) y(2,5) en una unidad y disminuye el flujo en los arcos (1,3) y(3,5) en una unidad. El arco (2,5) es el arco de bloqueo único y, por lo tanto, se selecciona para salir del árbol. Al soltar el arco (2,5) del árbol se obtienen dos subárboles: 𝑇1que consta de los nodos 1, 2, 3, 4 y 𝑇2que consta de los nodos 5 y 6. Introduciendo el arco (3,5), se obtiene nuevamente un árbol de expansión, véase la Figura 3.6 (c). En este árbol de expansión, los potenciales de nodo de los nodos 5 y 6 son una unidad menor que los del árbol de expansión anterior. En la solución de árbol de expansión factible que se muestra en la Figura 3.6 (c), 𝐿= {(2,3),(5,4)} y𝑈= {(2,5),(4,6)}. En esta solución, el arco (4,6) es el único arco elegible: su violación es igual a una unidad. Por lo tanto, se introduce el arco (4,6) en el árbol. La Figura 3.6 (c) muestra el ciclo resultante y su orientación. Puede aumentarse una unidad de flujo adicional a lo largo de la orientación de este ciclo. Al enviar este flujo, se encuentra que el arco (3,5) es un arco de bloqueo, por lo que se elimina este arco del árbol de expansión actual. La Figura 3.6 (d) muestra el nuevo árbol de expansión. Como puede verificarse, esta solución no tiene un arco elegible y, por lo tanto, el algoritmo simplex en redes termina con esta solución. 3. problema de flujo de mínimo coste: algoritmo simplex en redes 65 (a) (c) (b) (d) Figura 3.6: Ejemplo numérico para el algoritmo simplex en redes. 3.7 Árbol de expansión altamente factible El algoritmo simplex en redes no necesariamente termina en un número finito de iteraciones a menos que se imponga alguna restricción adicional sobre la elección de los arcos de entrada y salida. Al mantener un tipo especial de árbol de expansión, el árbol de expansión altamente factible, el algoritmo simplex en redes termina de manera finita; además, en la práctica se ejecuta más rápido. Sea (𝑇 , 𝐿, 𝑈)una estructura de árbol de expansión para un problema de flujo de mínimo coste con datos integrales. Los arcos de los árboles apuntan hacia arriba (hacia la raíz) o apuntan hacia abajo (dirección contraria a la raíz). Se establecen dos defini- 72 trabajo fin de máster Figura 3.10: Cálculo de los potenciales de nodo. ruta más corta no tienen cotas superiores, cada arco que no pertenece al árbol tiene el valor de su cota inferior. El algoritmo selecciona un arco (𝑘, 𝑙)no perteneciente al árbol con un coste reducido negativo para introducirlo en el árbol de expansión. La adición del arco (𝑘, 𝑙)al árbol crea un ciclo que se orienta en la misma dirección que el arco (𝑘, 𝑙). Sea 𝑤el vértice esquina de este ciclo. La Figura 3.11 ilustra esto. En este ciclo, cada arco desde el nodo 𝑙al nodo 𝑤es un arco hacia atrás y cada arco desde el nodo 𝑤al nodo 𝑘es un arco hacia adelante. En consecuencia, el arco saliente estaría en el segmento del nodo 𝑙al nodo 𝑤. De hecho, el arco saliente sería el arco (𝑝𝑟𝑒𝑑(𝑙), 𝑙) porque este arco tiene el valor de flujo más pequeño entre todos los arcos en el segmento desde el nodo 𝑙al nodo 𝑤. Luego, el algoritmo aumentaría los potenciales de los nodos en el subárbol partiendo del nodo 𝑙con una cantidad de |𝑐𝜋 𝑘𝑙|, actualizaría los índices del árbol y repetiría los cálculos hasta que todos los arcos ajenos al árbol tengan costes reducidos no negativos. Cuando el algoritmo termina, el árbol final sería un árbol de ruta más corta (un árbol en el que la ruta dirigida desde los nodos a todos los demás nodos es la ruta más corta). Recuérdese, que al implementar el algoritmo simplex en redes para el problema de flujo de mínimo coste, se mantiene los valores de flujo para todos los arcos porque se necesitan estos valores para identificar el arco saliente. Sin embargo, para el problema de la ruta más corta, puede determinarse el arco saliente sin considerar los valores de flujo. Si (𝑘, 𝑙)es el arco de entrada, entonces (𝑝𝑟𝑒𝑑(𝑙), 𝑙)es el arco de salida. Por lo tanto, el algoritmo simplex en redes para el problema de la ruta más corta no necesita 3. problema de flujo de mínimo coste: algoritmo simplex en redes 73 Figura 3.11: Elección del arco saliente. mantener los flujos de arco. Además, la actualización de los índices del árbol es más sencilla para el problema de la ruta más corta. El algoritmo simplex en redes aplicado al problema del camino más corto es similar al Algoritmo 6. Sea 𝑑(𝑖)=−𝜋(𝑖)la etiqueta de distancia válida, que representa la distancia de un camino dirigido desde la fuente 𝑠hasta el nodo 𝑖. El algoritmo simplex en redes busca un arco que satisfaga la condición 𝜋(𝑗)< 𝜋(𝑖) − 𝑐𝑖𝑗, y establece la etiqueta de distancia válida del nodo 𝑗igual a 𝑑(𝑖)+𝑐𝑖𝑗. En cada iteración, el algoritmo simplex en redes selecciona un arco (𝑖, 𝑗)con 𝑐𝜋 𝑖𝑗 <0, siendo 𝑐𝜋 𝑖𝑗 =𝑐𝑖𝑗 −𝜋(𝑖) + 𝜋(𝑗) = 𝑐𝑖𝑗 +𝑑(𝑖)−𝑑(𝑗). Luego, el algoritmo aumenta el potencial de cada nodo en el subárbol que parte del nodo 𝑗en una cantidad |𝑐𝜋 𝑖𝑗|que equivale a disminuir la etiqueta de distancia válida de todos los nodos en el subárbol enraizado en el nodo 𝑗en una cantidad |𝑐𝜋 𝑖𝑗|. 74 trabajo fin de máster 𝑑(𝑠)←0y𝑝𝑟𝑒𝑑(𝑠)←0; 𝑑(𝑗)←∞para cada 𝑗∈𝑁− {𝑠}; mientras algún arco (i, j) satisfaga d(j) > d(i) + c𝑖𝑗 hacer 𝑑(𝑗)←𝑑(𝑖) + 𝑐𝑖𝑗; 𝑝𝑟𝑒𝑑(𝑗)←𝑖; fin Algoritmo 6: Algoritmo de corrección del etiquetado. Si la red no contiene ningún ciclo negativo, el algoritmo simplex en redes termina con un árbol de ruta más corto. Cuando la red contiene algún ciclo negativo, el algoritmo eventualmente encuentra una situación como la que se muestra en la Figura 3.12. Este tipo de situación ocurre sólo cuando el nodo inicial del arco de entrada (𝑘, 𝑙) pertenece a 𝐷(𝑙), el conjunto de descendientes del nodo 𝑙. El algoritmo simplex en redes puede detectar esta situación fácilmente sin ningún aumento significativo en su esfuerzo computacional: después de introducir un arco (𝑘, 𝑙), el algoritmo actualiza los potenciales de todos los nodos en 𝐷(𝑙); en ese momento, puede verificar si 𝑘∈𝐷(𝑙)y, de ser así, terminar. En este caso, rastrear los índices predecesores producirá un ciclo negativo. La versión genérica del algoritmo simplex en redes para el problema de la ruta más corta se ejecuta en un tiempo pseudopolinomial. Este resultado se deriva de que: 1. para cada nodo 𝑖,−𝑛𝐶 ≤𝜋(𝑖)≤𝑛𝐶 (porque la longitud de cada camino dirigido desde 𝑠a el nodo 𝑖se encuentra entre −𝑛𝐶 y𝑛𝐶), y 2. cada iteración aumenta el valor de al menos un potencial de nodo. Hay, sin embargo, implementaciones que resuelven en tiempo polinomial, alguna de ellas son la regla de pivote del primer arco elegible, la regla de pivote de Dantzig o la regla del pivote escalado. 3.7.2 Algoritmo simplex en redes para el problema de flujo máximo El problema de flujo máximo puede considerarse una versión particular del problema de flujo de mínimo coste, obtenido al introducir un arco adicional (𝑡, 𝑠)con coeficiente de coste 𝑐𝑡𝑠 = −1 y una cota superior 𝑢𝑡𝑠 = ∞, y estableciendo 𝑐𝑖𝑗 = 0 para todos los arcos (𝑖, 𝑗)en 𝐴. Para simplificar la notación, de ahora en adelante asúmase que 𝐴 3. problema de flujo de mínimo coste: algoritmo simplex en redes 75 Figura 3.12: Detectando un ciclo negativo en la red. representa el conjunto 𝐴∪ {(𝑡, 𝑠)}. La formulación del problema es pues: minimizar −𝑥𝑡𝑠,(3.30) s. a ∑ {𝑗∶(𝑖,𝑗)∈𝐴} 𝑥𝑖𝑗 −∑ {𝑗∶(𝑗,𝑖)∈𝐴} 𝑥𝑗𝑖 = 0 para todo 𝑖∈𝑁, (3.31) 0≤𝑥𝑖𝑗 ≤𝑢𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴. (3.32) Minimizar −𝑥𝑡𝑠 equivale a maximizar 𝑥𝑡𝑠, lo que equivale a maximizar el flujo neto enviado desde la fuente al sumidero, ya que este flujo regresa a la fuente a través del arco (𝑡, 𝑠). Esto explica por qué la entrada es igual a la salida para todos los nodos de la red, incluidos los nodos fuentes y sumideros. En cualquier solución de árbol de expansión factible que lleve una cantidad positiva de flujo desde la fuente al sumidero (𝑥𝑡𝑠 >0), el arco (𝑡, 𝑠)debe de estar en el árbol de expansión. En consecuencia, el árbol de expansión para el problema de flujo máximo consta de dos subárboles de 𝐺unidos por el arco (𝑡, 𝑠)(ver Figura 3.13). Sean 𝑇𝑠y𝑇𝑡 los subárboles que contienen los nodo 𝑠y𝑡. 76 trabajo fin de máster Figura 3.13: Árbol de expansión para el problema de máximo flujo. Los potenciales de nodo correspondientes a un árbol de expansión factible del problema de flujo máximo se obtienen de la siguiente manera. Dado que puede establecerse el potencial de nodo arbitrariamente, sea 𝜋(𝑡)=0. Además, puesto que el coste reducido del arco (𝑡, 𝑠)debe ser cero, 0 = 𝑐𝜋 𝑡𝑠 =𝑐𝑡𝑠 −𝜋𝑡+𝜋𝑠= −1 + 𝜋𝑠, luego, 𝜋𝑠= 1. Dado que 1. el coste reducido de cada arco en 𝑇𝑠y𝑇𝑡debe ser cero, y 2. los costes de estos arcos son también cero, los potenciales de nodo tienen los siguientes valores: 𝜋(𝑖)=1para todo nodo 𝑖∈𝑇𝑠, y𝜋(𝑖)=0para todo nodo 𝑖∈𝑇𝑡. 3. problema de flujo de mínimo coste: algoritmo simplex en redes 77 Cada solución de árbol de expansión del problema de flujo máximo define un corte 𝑠−𝑡 ([𝑆, 𝑆]) en la red original obtenida estableciendo 𝑆=𝑇𝑠y𝑆=𝑇𝑡. Cada arco de este corte es un arco sin árbol; su flujo tiene valor cero o es igual a la capacidad del arco. Para cada arco hacia adelante (𝑖, 𝑗)en el corte, 𝑐𝜋 𝑖𝑗 = −1, y para cada arco hacia atrás (𝑖, 𝑗)en el corte, 𝑐𝜋 𝑖𝑗 = 1. Además, para todo arco (𝑖, −𝑗)que no está en el corte, 𝑐𝜋 𝑖𝑗 = 0. En consecuencia, si cada arco hacia adelante en el corte tiene un valor de flujo igual a la capacidad del arco y cada arco hacia atrás tiene flujo cero, esta solución de árbol de expansión satisface las condiciones de optimalidad (3.18)-(3.20) y, por lo tanto, debe ser óptima. Por otro lado, si en la solución actual del árbol de expansión, algún arco hacia adelante en el corte tiene un flujo nulo o el flujo en algún arco hacia atrás es igual a la capacidad del arco, todos estos arcos tienen una violación unitaria. Por lo tanto, puede seleccionarse cualquiera de estos arcos para ingresar al árbol de expansión. Supóngase que se selecciona el arco (𝑘, 𝑙). La introducción de este arco en el árbol crea un ciclo que contiene el arco (𝑡, 𝑠)como un arco hacia adelante (vea la Figura 3.14). El algoritmo aumenta el flujo máximo posible en este ciclo e identifica un arco de bloqueo. Dejar caer este arco nuevamente crea dos subárboles unidos por el arco (𝑡, 𝑠). Este nuevo árbol constituye un árbol de expansión para la próxima iteración. Este algoritmo es un algoritmo de camino incremental: la estructura de árbol permite determinar la ruta desde la fuente hasta el sumidero con mucha facilidad. En este sentido, el algoritmo simplex en redes tiene una ventaja sobre otros tipos de algoritmos de camino incremental para el problema de flujo máximo. Sin embargo, debido a la degeneración, es posible que el algoritmo simplex en redes no envíe una cantidad positiva de flujo desde la fuente al sumidero en cada iteración. Teorema 3.6 (Teorema de flujo máximo y mínimo corte). El valor máximo del flujo desde la fuente 𝑠hasta el sumidero 𝑡en una red capacitada es igual a la mínima capacidad entre todos los cortes 𝑠−𝑡. El algoritmo simplex en redes para el problema de flujo máximo proporciona una prueba del Teorema 3.6. El algoritmo termina cuando se capacita cada arco hacia adelante en el corte y cada arco hacia atrás tiene un flujo nulo. Esta condición de terminación implica que el valor de flujo máximo actual es igual a la capacidad del corte 𝑠−𝑡definido por los subárboles 𝑇𝑠y𝑇𝑡, y por lo tanto el valor de un flujo máximo desde el nodo 𝑠al nodo 𝑡es igual a la capacidad del corte mínimo 𝑠−𝑡. Así como la mecánica del algoritmo simplex en redes se vuelve más simple en el contexto del problema de flujo máximo, también lo hace en el contexto de un árbol de expansión altamente factible. Si se designa el sumidero como nodo raíz, la definición de un árbol de expansión altamente factible implica que pueda enviarse una cantidad 78 trabajo fin de máster Figura 3.14: Formación de un ciclo. positiva de flujo desde cada nodo en 𝑇𝑡al nodo sumidero 𝑡sin violar ninguno de las restricciones de flujo. Por lo tanto, todo arco en 𝑇𝑡cuyo valor de flujo es cero debe apuntar hacia el nodo sumidero 𝑡y todo arco en 𝑇𝑡cuyo valor de flujo, que es igual a la cota superior del arco, debe apuntar en sentido contrario al nodo 𝑡. Además, el criterio del arco saliente se reduce a seleccionar un arco de bloqueo en el ciclo de pivote que está más alejado del nodo sumidero cuando se atraviesa el ciclo en la dirección del arco (𝑡, 𝑠)comenzando desde el nodo 𝑡. Cada pivote degenerado selecciona un arco entrante que incide en algún nodo de 𝑇𝑡. La observación anterior implica que cada arco de bloqueo debe ser un arco en 𝑇𝑠. En consecuencia, cada pivote degenerado aumenta el tamaño de 𝑇𝑡, por lo que el algoritmo puede realizar como máximo 𝑛pivotes degenerados consecutivos. Podría darse que el problema de flujo de coste mínimo no satisficiera esto: para el problema más general, el número de pivotes degenerados consecutivos puede ser exponencialmente grande. La discusión anterior muestra que cuando se implementa para mantener un árbol de 3. problema de flujo de mínimo coste: algoritmo simplex en redes 79 expansión altamente factible, el algoritmo simplex en redes realiza (𝑛2𝑈)iteraciones para el problema de flujo máximo. Este resultado se deriva del hecho de que el número de pivotes no degenerados es como máximo 𝑛𝑈, cota superior del valor de flujo máximo. Esta magnitud en el número de iteraciones no es polinomial, por lo que no es satisfactorio desde una perspectiva teórica. 3.7.3 Análisis de sensibilidad El propósito del análisis de sensibilidad es determinar cambios en la solución óptima del problema de flujo de mínimo coste como producto de cambios en los datos (vector de oferta/demanda, capacidad o coste de cualquier arco). El análisis de sensibilidad adopta el siguiente enfoque básico. Primero se determina el efecto de un cambio dado en los datos sobre la viabilidad y la optimización de la solución asumiendo que la estructura del árbol de expansión permanece sin cambios. Si el cambio afecta a la optimalidad de la estructura del árbol de expansión, se realizan pivotes (primarios) para lograr la optimización. Siempre que el cambio destruye la viabilidad de la estructura del árbol de expansión, se realizan pivotes duales para lograr su viabilidad. Sea 𝑥∗una solución óptima del problema de flujo de mínimo coste, (𝑇∗, 𝐿∗, 𝑈∗)su correspondiente estructura de árbol de expansión y 𝜋∗denota los potenciales de nodo correspondientes. Primero se considera el análisis de sensibilidad con respecto a los cambios en los coeficientes de coste. Coste del análisis de sensibilidad Supóngase que el coste de un arco (𝑝, 𝑞)aumenta en 𝜆unidades. El análisis será distinto si el arco (𝑝, 𝑞)está en el árbol o no. Caso 1. El arco (𝑝, 𝑞)no está en el árbol. En este caso, cambiar el coste del arco (𝑝, 𝑞)no cambia los potenciales de nodo de la estructura del árbol de expansión. El coste reducido modificado del arco (𝑝, 𝑞)es 𝑐𝜋∗ 𝑝𝑞 +𝜆. Si el coste reducido modificado satisface la condición (3.19) o (3.20) (la que sea apropiada) la estructura del árbol de expansión sigue siendo óptima. De lo contrario, vuelve a optimizarse la solución utilizando el algoritmo simplex en redes con (𝑇∗, 𝐿∗, 𝑈∗)como estructura de árbol de expansión inicial. 80 trabajo fin de máster Caso 2. El arco (𝑝, 𝑞)está en el árbol. En este caso, cambiar el coste del arco (𝑝, 𝑞)cambia algunos potenciales de nodo. Si el arco (𝑝, 𝑞)es un arco que apunta hacia arriba en el árbol de expansión, los potenciales de todos los nodos en 𝐷(𝑝)aumentan en 𝜆, y si (𝑝, 𝑞)es un arco que apunta hacia abajo, los potenciales de todos los nodos en 𝐷(𝑞)disminuyen en 𝜆. Nótese que estos cambios alteran los costes reducidos de aquellos arcos ajenos al árbol que pertenecen al corte ([𝐷(𝑞), 𝐷(𝑞)]). Si todos los arcos que no son de árbol aún satisfacen la condición de optimalidad, la estructura actual del árbol de expansión sigue siendo óptima; de lo contrario, vuelve a optimizarse la solución utilizando el algoritmo simplex en redes. Análisis de sensibilidad de oferta/demanda Supóngase que el vector de oferta/demanda 𝑏(𝑘)del nodo 𝑘aumenta en 𝜆y la oferta/demanda 𝑏(𝑙)de otro nodo 𝑙disminuye en 𝜆(dado que: ∑𝑖∈𝑁𝑏(𝑖)=0, los suministros de dos nodos deben cambiar simultáneamente, en igual magnitud y direcciones opuestas). Las restricciones de balance de flujo requieren que deba enviarse 𝜆unidades de flujo desde el nodo 𝑘al nodo 𝑙. Sea 𝑃la ruta de árbol única del nodo 𝑘al nodo 𝑙. Sean 𝑃y𝑃, los conjuntos de arcos en 𝑃con la misma y opuesta dirección a la trayectoria respectivamente. El cambio de flujo máximo 𝛿𝑖𝑗 en un arco (𝑖, 𝑗) ∈ 𝑃que respeta las restricciones de flujo es 𝛿𝑖𝑗 ={𝑢𝑖𝑗 −𝑥𝑖𝑗 si (𝑖, 𝑗) ∈ 𝑃 , 𝑥𝑖𝑗 si (𝑖, 𝑗) ∈ 𝑃 . (3.33) Sea 𝛿= m ın{𝛿𝑖𝑗 ∶ (𝑖, 𝑗) ∈ 𝑃}.(3.34) Si 𝜆≤𝛿, se envía 𝜆unidades de flujo desde el nodo 𝑘al nodo 𝑙a lo largo de la ruta 𝑃. La solución modificada es factible para el problema modificado y dado que la modificación en 𝑏(𝑖)no afecta a la optimalidad de la solución, la solución resultante debe ser una solución óptima del problema modificado. Si 𝜆 > 𝛿, no puede enviarse 𝜆unidades de flujo desde el nodo 𝑘al nodo 𝑙a lo largo de los arcos del árbol de expansión y preservar su viabilidad. En este caso se envía 𝛿 unidades de flujo a lo largo de 𝑃y se reduce 𝜆a𝜆−𝛿. Sea 𝑥′el flujo actualizado, se realiza un pivote dual para obtener un nuevo árbol de expansión que podría permitir 3. problema de flujo de mínimo coste: algoritmo simplex en redes 81 el envío de flujo adicional desde el nodo 𝑘al nodo 𝑙a lo largo de la ruta del árbol. En un pivote dual, primero se decide la variable saliente y luego se identifica una variable entrante. Sea (𝑝, 𝑞)un arco en 𝑃que impide enviar flujo adicional desde el nodo 𝑘al nodo 𝑙. Si (𝑝, 𝑞) ∈ 𝑃, entonces 𝑥′ 𝑝𝑞 =𝑢𝑝𝑞 y si (𝑝, 𝑞) ∈ 𝑃, entonces 𝑥′ 𝑝𝑞 = 0. Se elimina el arco (𝑝, 𝑞)del árbol de expansión dividiendo el conjunto de nodos en dos subárboles. Sea 𝑆el subárbol que contiene el nodo 𝑘y𝑆el subárbol que contiene el nodo 𝑙. Dado que se desea enviar flujo adicional a través del corte [𝑆, 𝑆], los arcos elegibles para entrar en el árbol serían los arcos hacia adelante en el corte en su cota inferior o los arcos hacia atrás en sus cortes superiores. Si la red no contiene arcos capacitados no puede enviarse ningún flujo adicional desde el nodo 𝑘al nodo 𝑙y el problema modificado es inviable. Si la red contiene arcos capacitados, entonces, entre estos arcos, se selecciona un arco, por ejemplo (𝑔, ℎ), cuyo coste reducido tiene la menor magnitud. Se introduce el arco (𝑔, ℎ)en el árbol de expansión y se actualizan los potenciales de nodo. Luego se intenta nuevamente enviar 𝜆′=𝜆−𝛿unidades de flujo desde el nodo 𝑘al nodo 𝑙en la ruta del árbol. Si se logra, se termina; de lo contrario, se envía el máximo flujo posible y se realiza otro doble pivote para obtener una nueva estructura de árbol de expansión. Se repiten estos cálculos hasta que se establece un flujo factible en la red o se descubre que el problema modificado no es factible. Análisis de sensibilidad de la capacidad Finalmente, se considera el análisis de sensibilidad con respecto a las capacidades del arco. Considérese el caso en el que la capacidad de un arco (𝑝, 𝑞)aumenta en 𝜆 unidades. Siempre que se aumenta la capacidad de cualquier arco, la solución óptima anterior sigue siendo factible; para determinar si esta solución sigue siendo óptima, deben verificar las condiciones de optimalidad (3.18) - (3.20). Esté el arco en el árbol o no lo esté, aumentar 𝑢𝑝𝑞 en 𝜆unidades no afecta la condición de optimalidad. Sin embargo, si el arco (𝑝, 𝑞)con el valor de su cota superior no pertenece al árbol y su capacidad aumenta en 𝜆unidades, la condición de optimalidad (3.20) dicta que se debe aumentar el flujo en el arco en 𝜆unidades. Hacerlo crea un exceso de unidades 𝜆en el nodo 𝑞y un déficit de unidades 𝜆en el nodo 𝑝. Para lograr la viabilidad, debe enviarse unidades 𝜆del nodo 𝑞al nodo 𝑝. Se logra este objetivo utilizando el método descrito anteriormente en el análisis de sensibilidad de oferta/demanda. 88 trabajo fin de máster Si se considera un servicio fijo 𝑘, este coste reducido es similar al coste reducido utilizado anteriormente (2.3) para el problema de flujo de mínimo coste (Capítulo 3) ; la diferencia es que ahora se suma el precio de la arista 𝑤𝑖𝑗 al coste 𝑐𝑘 𝑖𝑗. Nótese como las restricciones de capacidad (4.2) relacionan las variables 𝑥𝑘 𝑖𝑗, sin ellas serían independientes. Los costes de las aristas 𝑤𝑖𝑗 relacionan los costes reducidos por servicio y evitan que sean independientes. Para caracterizar las soluciones óptimas del problema del flujo multiservicio, primero se escribe el problema dual del flujo multiservicio (4.1)-(4.4): Maximizar −∑ (𝑖,𝑗)∈𝐴 𝑢𝑖𝑗𝑤𝑖𝑗 + 𝑘 ∑ 𝑘=1 𝑏𝑘𝜋𝑘(4.7) s. a 𝑐𝜋,𝑘 𝑖𝑗 =𝑐𝑘 𝑖𝑗 +𝑤𝑖𝑗 −𝜋𝑘(𝑖) + 𝜋𝑘(𝑗)≥0para todo (𝑖, 𝑗) ∈ 𝐴y todo 𝑘= 1,…, 𝐾, (4.8) 𝑤𝑖𝑗 ≥0para todo (𝑖, 𝑗) ∈ 𝐴. (4.9) Las condiciones de optimalidad de la programación lineal, llamadas condiciones (óptimas) de holgura complementaria, establecen que una solución factible primal 𝑥y una solución factible dual (𝑤, 𝜋𝑘)son óptimas para los problemas respectivos si y sólo si el servicio de cada variable primal (dual) y la holgura en la correspondiente restricción dual (primal) es cero. Las condiciones de holgura complementaria para el par primaldual del problema de flujo multiservicio asumen la siguiente forma especial. (Aquí se usa 𝑦𝑘 𝑖𝑗 para denotar un valor específico de la variable de flujo 𝑥𝑘 𝑖𝑗.) Teorema 4.1 (Condiciones de holgura del problema de flujo multiservicio). Los flujos de cada servicio 𝑦𝑘 𝑖𝑗 son óptimos en el problema de flujo multiservicio (4.1)- (4.4) con cada 𝑢𝑘 𝑖𝑗 = +∞ si y sólo si son factibles para alguna elección de costes (no negativos) de aristas 𝑤𝑖𝑗 y potenciales de nodos (sin restricciones en el signo) 𝜋𝑘(𝑖), los costes reducidos y los flujos de las aristas satisfacen las siguientes condiciones de holgura: (a) 𝑤𝑖𝑗 (∑ 1≤𝑘≤𝐾 𝑦𝑘 𝑖𝑗 −𝑢𝑖𝑗)= 0 para toda arista (𝑖, 𝑗) ∈ 𝐴. (4.10) (b) 𝑐𝜋,𝑘 𝑖𝑗 ≥0para toda arista (𝑖, 𝑗) ∈ 𝐴y(4.11) todos los servicios 𝑘= 1,2,…, 𝐾. (4.12) (c) 𝑐𝜋,𝑘 𝑖𝑗 𝑦𝑘 𝑖𝑗 = 0 para toda arista (𝑖, 𝑗) ∈ 𝐴y(4.13) todos los servicios 𝑘= 1,2,…, 𝐾. (4.14) 4. problema de flujo multiservicio 89 Los costes de aristas óptimos y potenciales de nodos óptimos son cualquier conjunto de costes de aristas y potenciales de nodos que satisfacen las condiciones de holgura. El siguiente teorema muestra la conexión entre los problemas de flujo multiservicio y monoservicio. Teorema 4.2 (Dualización parcial). Sea 𝑦𝑘 𝑖𝑗 los flujos óptimos y 𝑤𝑖𝑗 los costes óptimos de las aristas para el problema del flujo de multiservicio (4.1)-(4.4). Luego, para cada servicio 𝑘, las variables de flujo 𝑦𝑘 𝑖𝑗 para (𝑖, 𝑗) ∈ 𝐴resuelven el siguiente problema de flujo de coste mínimo (no capacitado): m ın {∑ (𝑖,𝑗)∈𝐴(𝑐𝑘 𝑖𝑗 +𝑤𝑖𝑗)𝑥𝑘 𝑖𝑗 ∶𝑥𝑘=𝑏, 𝑥𝑘 𝑖𝑗 ≥0para todo (𝑖, 𝑗) ∈ 𝐴}.(17.3) 4.3 Relajación Lagrangiana Para aplicar la relajación lagrangiana al problema de flujo multiservicio, se asocian multiplicadores 𝑤𝑖𝑗 lagrangianos no negativos a las restricciones (4.2), dando lugar al siguiente subproblema lagrangiano: 𝐿(𝑤)=m ın ∑ 1≤𝑘≤𝐾 𝑐𝑘𝑥𝑘+∑ (𝑖,𝑗)∈𝐴 𝑤𝑖𝑗 (∑ 1≤𝑘≤𝐾 𝑥𝑘 𝑖𝑗 −𝑢𝑖𝑗),(4.15) o, de forma equivalente, 𝐿(𝑤)=m ın ∑ 1≤𝑘≤𝐾∑ (𝑖,𝑗)∈𝐴(𝑐𝑘 𝑖𝑗 +𝑤𝑖𝑗)𝑥𝑘 𝑖𝑗 −∑ (𝑖,𝑗)∈𝐴 𝑤𝑖𝑗𝑢𝑖𝑗 (4.16) s. a 𝑥𝑘=𝑏𝑘para todo 𝑘= 1,…, 𝐾, (4.17) 𝑥𝑘 𝑖𝑗 ≥0para todo (𝑖, 𝑗) ∈ 𝐴y todo 𝑘= 1,2,…, 𝐾. (4.18) Nótese que dado el término −∑(𝑖,𝑗)∈𝐴𝑤𝑖𝑗𝑢𝑖𝑗 en la función objetivo del subproblema lagrangiano es siempre una constante. Sean cuales fuere los multiplicadores lagrangianos, este término es constante y, por tanto, se puede ignorar. En la función objetivo resultante del subproblema lagrangiano cada variable 𝑥𝑘 𝑖𝑗 tiene un coste 𝑐𝑘 𝑖𝑗 +𝑤𝑖𝑗. Dado que en este problema ninguna de las restricciones contienen las variables de flujo de más de un servicio, el problema se puede descomponer en problemas de flujo de mínimo coste separados entre sí, uno para cada servicio. En consecuencia, aplicando el procedimiento de optimización de subgradiante de la Sección 2.2.6, se puede alternativamente: 90 trabajo fin de máster 1. Hallar el conjunto de flujos de mínimo coste (para valores fijos de los multiplicadores lagrangianos 𝑤) con los coeficientes de coste 𝑐𝑘 𝑖𝑗 +𝑤𝑖𝑗. 2. Actualizar los multiplicadores con el algoritmo descrito en la Sección 2.2.6. En este último caso, si 𝑦𝑘 𝑖𝑗 denota la solución optima del subproblema de flujo de coste mínimo cuando los multiplicadores lagrangianos tienen el valor 𝑤𝑞 𝑖𝑗 en la 𝑞-ésima iteración, actualizando la fórmula del subgradiante queda: 𝑤𝑞+1 𝑖𝑗 =[𝑤𝑞 𝑖𝑗 +𝜃𝑞(∑ 1≤𝑘≤𝐾 𝑦𝑘 𝑖𝑗 −𝑢𝑖𝑗)]+ .(4.19) La notación [𝛼]+denota la parte positiva de 𝛼osea, m ax(𝛼, 0). El escalar 𝜃𝑞es el tamaño del paso que especifica como de lejos se mueve uno de la solución 𝑤𝑞 𝑖𝑗. Nótese que la formula actualizadora hace lo siguiente: 1. Si la solución 𝑦𝑘 𝑖𝑗 del subproblema usa más de la capacidad 𝑢𝑖𝑗 de la arista, incrementa el multiplicador 𝑤𝑞 𝑖𝑗 una cantidad (∑1≤𝑘≤𝐾𝑦𝑘 𝑖𝑗 −𝑢𝑖𝑗)en la arista (𝑖, 𝑗). 2. Si la solución 𝑦𝑘 𝑖𝑗 del subproblema usa menos de la capacidad disponible en la arista, reduce el multiplicador de lagrangiano una cantidad (∑1≤𝑘≤𝐾𝑦𝑘 𝑖𝑗)en la arista (𝑖, 𝑗). 3. Sin embargo, si este decrecimiento hace que el multiplicador 𝑤𝑞+1 𝑖𝑗 se haga negativo, entonces se reduce el valor del multiplicador sólo a cero. Se elige el tamaño del paso 𝜃𝑞para las interaciones 𝑞= 1,2,…, 𝐾, de acuerdo con el procedimiento descrito en la Sección 2.2.6. El método del subgradiante (Sección 2.2.6) para el problema de flujo multiservicio es un procedimiento para resolver el problema lineal (4.1)-(4.4) que es capaz de hacer uso de la relajación de las restricciones de balance de flujo. En el Teorema 2.2 se señaló que siempre que se aplica la relajación lagrangiana a cualquier problema lineal, como es el caso del problema de flujo multiservicio, el valor óptimo 𝐿∗= m ax𝑤≥0𝐿(𝑤) del problema dual lagrangiano es igual al valor óptimo de la función objetivo 𝑧∗del problema lineal. El uso de la optimización subgradiante para resolver el problema de multiplicadores lagrangianos es atractivo por varias razones. Primero, como se acaba de señalar, este planteamiento permite explotar la estructura de la red de flujo. En segundo lugar, las 4. problema de flujo multiservicio 91 fórmulas para actualizar los multiplicadores lagrangianos 𝑤𝑖𝑗 computacionalmente son sencillas y muy fáciles de codificar en un programa de ordenador. Sin embargo, este tipo de solución también tiene algunas limitaciones. Para asegurar la convergencia, son necesarios pasos pequeños; como resultado, el método no converge muy rápido. En segundo lugar, el método es de base dual y, por lo tanto, aunque converge a las variables duales óptimas 𝑤𝑖𝑗, las soluciones óptimas 𝑦𝑘 𝑖𝑗 de los subproblemas no tienen por qué coincidir con la solución óptima del problema de flujo multiservicio. De hecho, aunque se ha demostrado que si se establecen los valores óptimos de los multiplicadores lagrangianos (es decir, las variables duales óptimas para las restricciones de capacidad en el problema lineal (4.1)-(4.4)), los flujos óptimos 𝑦𝑘 𝑖𝑗 resuelven el subproblema lagrangiano (4.15)-(4.18). Estos subproblemas también pueden tener otras soluciones óptimas que no satisfacen las restricciones de capacidad. Supóngase que se resuelve el subproblema de Lagrange para cualquier elección de los multiplicadores mediante el algoritmo simplex en redes (Capítulo 3). Al hacerlo, se obtiene una solución óptima de árbol generador para cada servicio, que corresponde a una base 𝑘de las restricciones (4.17) de capacidad de la red. Se puede extender esta base a una base general del problema lineal (4.1)-(4.4), como se muestra en la Figura 4.1. La matriz de identidad en la fila superior de esta matriz corresponde a las variables de holgura 𝑠𝑖𝑗 en la igualdad (4.5) de las restricciones de capacidad. Obsérvese que en la solución correspondiente a la base : Cada variable de flujo 𝑥𝑘para 𝑘= 1,2,…, 𝐾 es igual a la solución 𝑦𝑘del subproblema lagrangiano. los valores de las variables de holgura vienen dados por 𝑠𝑖𝑗 =𝑢𝑖𝑗 −∑1≤𝑘≤𝐾𝑦𝑘 𝑖𝑗. Esta base es factible si todas las variables de holgura son no negativas. Si alguna de ellas es negativa la base es inviable. En cualquier caso, se puede aplicar métodos de programación lineal estándar (ya sea el método simplex primal o el método simplex dual) usando esta base como una solución inicial para resolver el problema. Esta base proporciona una solución inicialmente atractiva para el método símplex porque permite explotar la estructura de red y los costes de flujo para generar una solución inicial ”inteligente”. Para problemas que están ”capacitados modestamente”, la base puede ser una muy buena aproximación de la solución óptima del problema y, por lo tanto, puede mejorar en gran medida el rendimiento del método simplex. 92 trabajo fin de máster coeficientes de la restricción de capacidad Servicio 1 Servicio 2 Servicio 3 coeficientes de las variables de flujo en las restricciones de capacidad Servicio K Figura 4.1: Base formado por las bases de los subproblemas (4.1)-(4.4) de cada servicio. A continuación, se considera un planteamiento alternativo, conocido como descomposición de Dantzig-Wolfe, para resolver el problema dual de Lagrange. Este algoritmo requiere considerablemente más trabajo (la solución de un problema lineal) en cada iteración para actualizar los multiplicadores de Lagrange, pero ha demostrado converger más rápido que el procedimiento de optimización subgradiente para varias clases de problemas. En lugar de describir el procedimiento de descomposición de Dantzig- Wolfe como una variación de la relajación lagrangiana, primero se desarrolla desde un punto de vista de programación lineal a gran escala que proporciona una perspectiva algo diferente. 4. problema de flujo multiservicio 93 4.4 Generación de columnas Considérese un caso especial del problema de flujo multiservicio: cada servicio 𝑘tiene un solo nodo fuente 𝑠𝑘y un solo nodo sumidero 𝑡𝑘y un flujo de 𝑑𝑘unidades. No se impone límite a la capacidad de flujo de los servicios individualmente salvo las restricciones de capacidad. Por lo tanto, para cada servicio 𝑘, las restricciones del subproblema 𝑥𝑘=𝑏𝑘,𝑥𝑘≥0definen un problema de camino más corto. En este modelo, para cualquier elección 𝑤𝑖𝑗 de los multiplicadores de Lagrange, la relajación lagrangiana de las restricciones de capacidad requiere la solución de un problema de ruta más corta por cada servicio. 4.4.1 Reformulación de flujos en ruta En la formulación de los problemas de flujo de red se puede definir el flujo en arcos o definir el flujo en caminos y en ciclos. Por ”flujo en arco” se entiende un vector 𝑥= {𝑥𝑖𝑗}que satisface las siguientes restricciones: ∑ {𝑗∶(𝑖,𝑗)∈𝐴} 𝑥𝑖𝑗 −∑ {𝑗∶(𝑗,𝑖)∈𝐴} 𝑥𝑗𝑖 = −𝑒(𝑖)para todo 𝑖∈𝑁, (4.20) 0≤𝑥𝑖𝑗 ≤𝑢𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴. (4.21) donde ∑𝑛 𝑖=1 𝑒(𝑖)=0. En este modelo se ha reemplazado la oferta/demanda 𝑏(𝑖)de cada nodo 𝑖por otro término, −𝑒(𝑖).𝑒(𝑖)es el desequilibrio del nodo y representa lo que entra menos lo que sale del nodo 𝑖. Si entra más de lo que sale 𝑒(𝑖)>0se dice que 𝑖es un nodo en exceso. Si sale más de lo que entra 𝑒(𝑖)<0se dice que 𝑖es un nodo en deficit. Si entra tanto como sale se dice que 𝑖es un nodo balanceado. Si 𝑒= −𝑏el flujo 𝑥es factible para el problema de flujo de mínimo coste. En la formulación de flujo de arco las variables de decisión básicas son los flujos 𝑥𝑖𝑗 en los arcos (𝑖, 𝑗) ∈ 𝐴. La formulación de flujos en camino o ciclo empieza con una enumeración de los caminos 𝑃entre cualquier par de nodos y la colección de todos los ciclos . La variable de decisión en la formulación del flujo en camino y del flujo en ciclo son 𝑓(𝑃)y𝑓(𝑊). Se define estas variables para cada camino 𝑝∈y cada ciclo 𝑊∈. El flujo 𝑥𝑖𝑗 en el arco (𝑖, 𝑗)es igual a la suma de los flujos 𝑓(𝑃)y𝑓(𝑊)para todos los caminos 𝑃y ciclos 𝑊que contienen a este arco. Sea 𝛿𝑖𝑗(𝑃)=1si (𝑖, 𝑗) ∈ 𝑃y 94 trabajo fin de máster 𝛿𝑖𝑗(𝑃)=0en caso contrario. Entonces 𝑥𝑖𝑗 =∑ 𝑃∈ 𝛿𝑖𝑗(𝑃)𝑓(𝑃) + ∑ 𝑊∈ 𝛿𝑖𝑗(𝑊)𝑓().(4.22) Por tanto, los flujos en caminos y los flujos en ciclos determinan de forma única los flujos en arcos. Teorema 4.3 (Teorema de la descomposición de flujo). Cada flujo en camino o en ciclo tiene una única representación como flujos no negativos en arcos. A su vez, con las siguientes dos propiedades cada flujo de arco no negativo 𝑥se puede representar como un flujo en camino y ciclo (no necesariamente de forma única): 1. Cada camino dirigido con flujo positivo conecta un nodo deficitario con un nodo en exceso. 2. Al menos 𝑛+𝑚caminos y ciclos tienen flujo no negativo; además, como mucho 𝑚 ciclos tienen flujo no negativo. Sea 𝑥un flujo para el cual 𝑒(𝑖)=0por cada 𝑖∈𝑁. Cuando se aplica el algoritmo de descomposición de flujo a una circulación, cada iteración describe un ciclo dirigido que consta únicamente de arcos con flujo positivo y, posteriormente, reduce el flujo en al menos un arco a cero. En consecuencia, una circulación se descompone en flujos a lo largo de como máximo 𝑚ciclos dirigidos. Propiedad 4.1.Una circulación 𝑥puede representarse como un flujo en ciclo a lo largo de la mayoría de los 𝑚ciclos dirigidos. El Teorema 4.3 tiene varias consecuencias importantes. Por ejemplo, permite comparar dos soluciones de un problema de flujo de red de una manera conveniente y mostrar cómo se puede construir una solución a partir de otra mediante una secuencia de operaciones simples. Para definir los ciclos de aumento es primero necesario introducir los arcos dirigidos hacia delante y los arcos dirigidos hacia atrás. Los arcos dirigidos hacia delante son aquellos que se corresponden con la dirección del flujo y los arcos dirigidos hacia atras aquellos contrarios a este, véase la Figura 4.2. Un ciclo 𝑊(no necesariamente dirigido) de 𝐺es un ciclo de aumento con respecto al flujo 𝑥si al aumentar una cantidad positiva de flujo 𝑓(𝑊)por el ciclo, el flujo sigue siendo factible. Este crecimiento aumenta el flujo en los arcos dirigidos hacia adelante en el ciclo 𝑊y disminuye el flujo en los arcos dirigidos hacia atrás en el 4. problema de flujo multiservicio 95 arcos hacia delante arcos hacia atrás W-W Figura 4.2: Ejemplos de arcos hacia delante y arcos hacia atrás en dos ciclos. ciclo. Por lo tanto, un ciclo 𝑊es un ciclo de aumento de 𝐺si 𝑥𝑖𝑗 ≤𝑢𝑖𝑗 para cada arco dirigido hacia adelante (𝑖, 𝑗)y𝑥𝑖𝑗 ≥0para cada arco dirigido hacia atrás (𝑖, 𝑗). A continuación, se amplía la notación de 𝛿𝑖𝑗(𝑊)para ciclos que no están necesariamente dirigidos. Se define 𝛿𝑖𝑗(𝑊)=1si (𝑖, 𝑗)es un arco dirigido hacia adelante en el ciclo 𝑊,𝛿𝑖𝑗(𝑊) = −1 si (𝑖, 𝑗)es un arco dirigido hacia atrás en el ciclo 𝑊, e igual a 0en caso contrario. En términos de redes residuales, cada ciclo de aumento 𝑊con respecto a un flujo 𝑥 corresponde a un ciclo dirigido 𝑊en 𝐺(𝑥), y viceversa. Se define el coste de un ciclo de aumento 𝑊como 𝑐(𝑊) = ∑(𝑖,𝑗)∈𝑊𝑐𝑖𝑗𝛿𝑖𝑗(𝑊). El coste de un ciclo de aumento representa el cambio en el coste de una solución factible si se aumenta en 1 unidad el flujo a lo largo del ciclo. El cambio en el coste de flujo para aumentar 𝑓(𝑊)unidades a lo largo del ciclo 𝑊es 𝑐(𝑊)𝑓(𝑊). Supóngase que 𝑥y𝑥0son dos soluciones factibles cualesquiera del problema de flujo con mínimo coste. Alguna circulación factible 𝑥1en 𝐺0fatisface la propiedad 𝑥= 𝑥0+𝑥1. La Propiedad 4.1 supone que se puede representar la circulación 𝑥1como flujos en ciclo 𝑓(𝑊1), 𝑓(𝑊2),…, 𝑓(𝑊𝑟), con 𝑟≤𝑚. Nótese que cada uno de los ciclos 𝑊1, 𝑊2,…, 𝑊𝑟es un ciclo de aumento en 𝐺(𝑥0). Además, ∑ (𝑖,𝑗)∈𝐴 𝑐𝑖𝑗𝑥𝑖𝑗 =∑ (𝑖,𝑗)∈𝐴 𝑐𝑖𝑗𝑥0 𝑖𝑗 +∑ (𝑖,𝑗)∈𝐺(𝑥0 𝑐𝑖𝑗𝑥1 𝑖𝑗 (4.23) =∑ (𝑖,𝑗)∈𝐴 𝑐𝑖𝑗𝑥0 𝑖𝑗 +∑ (𝑖,𝑗)∈𝐺(𝑥0 𝑐𝑖𝑗 [𝑟 ∑ 𝑘=1 𝛿𝑖𝑗(𝑊𝑘)𝑓(𝑊𝑘)](4.24) =∑ (𝑖,𝑗)∈𝐴 𝑐𝑖𝑗𝑥0 𝑖𝑗 + 𝑟 ∑ 𝑘=1 𝑐(𝑊𝑘)𝑓(𝑊𝑘).(4.25) Teorema 4.4 (Teorema del ciclo de aumento). Sean 𝑥y𝑥0dos soluciones posi- 96 trabajo fin de máster bles de un problema de flujo en red. Entonces 𝑥es igual a 𝑥0más el flujo en como máximo 𝑚ciclos dirigidos de 𝐺(𝑥0). Además, el coste de 𝑥es igual al coste de 𝑥0más el coste del flujo en estos ciclos de aumento. Teorema 4.5 (Teorema de Optimidad de Ciclo Negativo). Una solución factible 𝑥∗del problema de flujo de mínimo coste es una solución óptima si y sólo si la red residual 𝐺(𝑥∗)no contiene un ciclo dirigido de coste negativo. Se reformula el problema de flujo multiservicio utilizando flujos de trayectoria y flujos de ciclo en lugar de flujos de arco. Para simplificar más, supóngase que para cada servicio el coste de cada ciclo 𝑊en la red subyacente es no negativo. El problema satisface esta condición, por ejemplo, si los costes de flujo en cada arco son no negativos. Si se impone esta condición de coste de ciclo no negativo, entonces en alguna de las soluciones óptimas del problema, el flujo en cada ciclo es cero, por lo que se puede eliminar las variables de flujo de los ciclos. Así se podrá representar cualquier solución potencialmente óptima como la suma de los flujos en las rutas dirigidas. Para cada servicio 𝑘, sea 𝑃𝑘la colección de todas las rutas dirigidas desde la fuente 𝑠𝑘hasta el sumidero 𝑡𝑘en la red 𝐺= (𝑁, 𝐴). En la formulación de flujo de ruta, cada variable de decisión 𝑓(𝑃)es el flujo en alguna ruta 𝑃para el 𝑘-ésimo servicio, se define esta variable para cada ruta 𝑃en 𝑃𝑘. Sea 𝛿𝑖𝑗(𝑃)una variable indicadora. Indica la trayectoria en arco, es decir, 𝛿𝑖𝑗(𝑃)es igual a 1 si la trayectoria 𝑃contiene al arc(𝑖, 𝑗), y es 0 en caso contrario. El Teorema 4.3 establece que un flujo en arco 𝑥𝑘 𝑖𝑗 óptimo siempre se puede descomponer en flujos de trayectoria 𝑓(𝑃)de la siguiente manera: 𝑓(𝑃) = 𝑥𝑘 𝑖𝑗 ∑ 𝑃∈𝑃𝑘 𝛿𝑖𝑗(𝑃)𝑓(𝑃).(4.26) Sea 𝑐𝑘(𝑃) = ∑(𝑖,𝑗)∈𝐴𝑐𝑘 𝑖𝑗𝛿𝑖𝑗 (𝑃) = ∑(𝑖,𝑗)∈𝑃𝑐𝑘 𝑖𝑗 que denota el coste de la unidad de flujo en la trayectoria 𝑃∈𝑃𝐾con respecto al servicio 𝑘. Para cada servicio 𝑘, si se sustituyen las variables de flujo del arco en la función objetivo, se intercambia el orden de las sumas y se reordenan términos, se tiene que: ∑ (𝑖,𝑗)∈𝐴 𝑐𝑘 𝑖𝑗𝑥𝑘 𝑖𝑗 =∑ (𝑖,𝑗)∈𝐴 𝑐𝑘 𝑖𝑗 [∑ 𝑃∈𝑃𝑘 𝛿𝑖𝑗(𝑃)𝑓(𝑃)]=∑ 𝑃∈𝑃𝑘 𝑐𝑘(𝑃)𝑓(𝑃).(4.27) Esto muestra que se puede expresar el coste de cualquier solución como el coste de los flujos de arco o el coste de los flujos de ruta. 4. problema de flujo multiservicio 97 Al sustituir las variables de ruta en la formulación de flujo multiservicio, se obtiene la siguiente formulación de flujo de ruta: Minimizar ∑ 1≤𝑘≤𝐾∑ 𝑃∈𝑃𝑘 𝑐𝑘(𝑃)𝑓(𝑃)(4.28) s. a ∑ 1≤𝑘≤𝐾∑ 𝑃∈𝑃𝑘 𝛿𝑖𝑗(𝑃)𝑓(𝑃)≤𝑢𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴, (4.29) ∑ 𝑃∈𝑃𝑘 𝑓(𝑃) = 𝑑𝑘para todo 𝑘= 1,2,…, 𝐾, (4.30) 𝑓(𝑃)≤0para todo 𝑘= 1,2,…, 𝐾 y todo 𝑃∈𝑃𝑘.(4.31) Al formular este problema, se ha usado Teorema 4.3 que establece que se puede descomponer cualquier arco de flujo factible del sistema 𝑥𝑘=𝑏𝑘en un conjunto de flujos de trayectoria y flujos de ciclo de tal manera que los trayectos satisfacen la condición de balance de flujo (4.30). La formulación de flujo de ruta del problema de flujo multiservicio tiene una estructura de restricciones muy simple. El problema tiene una única restricción (4.29) para cada arco (𝑖, 𝑗)que establece que la suma de los flujos de la trayectoria que atraviesan el arco es como máximo 𝑢𝑖𝑗 (la capacidad del arco). Además, el problema tiene una única restricción (4.30) para cada servicio 𝑘que establece que el flujo total en todas las rutas que conectan la fuente 𝑠𝑘y el sumidero 𝑡𝑘del servicio 𝑘debe ser igual a la demanda 𝑑𝑘de este servicio. Para una red con 𝑛nodos, 𝑚arcos y 𝐾servicios, la formulación de flujo de ruta contiene 𝑚+𝐾restricciones (además de las restricciones de no negatividad impuestas a los valores de flujo de ruta (4.31)). Por el contrario, la formulación en arco (4.28) contiene 𝑚+𝑛𝐾 restricciones puesto que contiene una restricción de balance de flujo para cada combinación de nodo y servicio. Cómo sólo aparece una ruta en una de las restricciones (4.30), se puede aplicar una versión especifica del método simplex, conocido como método simplex de cota superior generalizada, para resolver la formulación de flujo de ruta de manera más eficiente. Este método esencialmente resuelve el problema como si sólo tuviera 𝑚restricciones de capacidad. Esta disminución en el número de restricciones tiene su precio, la formulación de flujo de ruta tiene una variable para cada ruta que, por servicio, conecta una fuente y un sumidero. El número de variables suele ser enorme, creciendo exponencialmente con el tamaño de la red. Pero se podría esperar que sólo muy pocos de los caminos recorran la solución óptima al problema. De hecho, se sabe que, como máximo, las 𝐾+𝑚trayectorias llevan un flujo positivo en alguna solución óptima al problema. Si se conociera el conjunto óptimo de caminos, o un muy buen conjunto de caminos, se 104 trabajo fin de máster como una mejora de la solución. El pricer resuelve el subproblema: m ın ∑ 1≤𝑘≤𝐾∑ 𝑃∈𝑃𝐾[𝑐𝑘(𝑃) + ∑ (𝑖,𝑗)∈𝑃 𝑤𝑖𝑗](4.42) s. a ∑ 𝑃∈𝑃𝑘 𝑓(𝑃) = 𝑑𝑘para todo 𝑘= 1,2,…, 𝐾, (4.43) 𝑓(𝑃)≥0para todo 𝑘= 1,2,…, 𝐾 y todo 𝑃∈𝑃𝐾.(4.44) donde 𝑐𝑘(𝑃) = ∑(𝑖,𝑗)∈𝑃𝑐𝑘 𝑖𝑗. Así cada pricer resuelve un problema de camino más corto con los pesos facilitados por el master y reportan si encuentran un camino más corto para su correspondiente servicio que el que está ya usando el master. El algoritmo mantiene una base de programación lineal e introduce una columna en cada iteracción. Descarta cualquier columna que abandona la base pero dado que generar columnas es un proceso costoso y estas columnas pueden tener un valor negativo en futuras iteraciones, se almacenan. Así, el algoritmo, en general, realiza más de un cambio en la base, resuelve el problema maestro restringido, halla los nuevos costes y se resuelven de nuevo los subproblemas de camino más corto. Cada subproblema de camino más corto se puede resolver de forma paralela. El algoritmo es finito, cada columna del problema maestro corresponde a un camino entre la cantidad finita de posibles caminos. Eventualmente el algoritmo generará todas las columnas del problema maestro. Los 𝐾subproblemas corresponden a la relajación lagrangiana con multiplicadores 𝑤𝑖𝑗 en cada arco (𝑖, 𝑗). El master fija los multiplicadores y resuelve el problema de multiplicadores de Lagrange. La descomposición de Dantzig-Wolfe es un método eficiente para resolver la relajación lagrangiana si se mide la eficiencia mediante el número de iteraciones de un algoritmo. El algoritmo de descomposición de Dantzig-Wolfe siempre mantiene una solución factible del problema. Como se vio en la Sección 4.5 la solución del subproblema aporta una cota de cómo de alejada está la solución factible de la óptima. Se puede parar el algoritmo en cualquier iteración y obtener no sólo una solución factible sino una medida de cuan alejado está la solución factible de la solución óptima. 4. problema de flujo multiservicio 105 4.7 Descomposición de la directiva de recursos El método de la directiva de recursos, en lugar de utilizar costes para descomponer el problema, asigna la capacidad de los arcos a los servicios. Cuando se aplica a la formulación del problema (17.1), asigna 𝑟𝑘 𝑖𝑗 ≤𝑢𝑘 𝑖𝑗 unidades de la capacidad 𝑢𝑖𝑗 del arco (𝑖, 𝑗)al servicio 𝑘, produciendo el siguiente problema de directiva de recursos: 𝑧= m ın ∑ 1≤𝑘≤𝐾 𝑐𝑘𝑥𝑘(4.45) ∑ 1≤𝑘≤𝐾 𝑟𝑘 𝑖𝑗 ≤𝑢𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴, (4.46) 𝑥𝑘=𝑏𝑘para 𝑘= 1,2,…, 𝐾, (4.47) 0≤𝑥𝑘 𝑖𝑗 ≤𝑟𝑘 𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴y todo 𝑘= 1,2,…, 𝐾. (4.48) La restricción (4.46) asegura que la asignación total de recursos para el arco (𝑖, 𝑗)no exceda la capacidad de ese arco. Sea 𝑟= (𝑟𝑘 𝑖𝑗)el vector de asignaciones de recursos. Propiedad 4.2.El problema de la directiva de recursos (4.45)-(4.48) es equivalente al problema original del flujo multiservicio (4.1)-(4.4) en el sentido de que 1. si (𝑥, 𝑟)es factible en el problema de la directiva de recursos, entonces 𝑥es factible y la función tiene el mismo valor objetivo que en el problema original, y 2. si 𝑥es factible en el problema original y 𝑟=𝑥, entonces (𝑥, 𝑟)es factible y la función objetivo tiene el mismo valor que en el problema de directiva de recursos. Véase el siguiente método para resolver el problema de la directiva de recursos (4.45)- (4.48). En lugar de resolver el problema eligiendo los vectores 𝑟y𝑥simultáneamente, se eligen secuencialmente. Primero se fijan las asignaciones de recursos 𝑟𝑘 𝑖𝑗 y luego se elige el flujo 𝑥𝑘 𝑖𝑗. Sea 𝑧(𝑟)el valor óptimo del problema de la directiva de recursos para un valor fijo de la asignación de recursos 𝑟y considérese el siguiente problema derivado de asignación de recursos: Minimizar 𝑧(𝑟)(4.49) s. a ∑ 1≤𝑘≤𝐾 𝑟𝑘 𝑖𝑗 ≤𝑢𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴, (4.50) 0≤𝑟𝑘 𝑖𝑗 ≤𝑢𝑘 𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴y todo 𝑘= 1,2,…, 𝐾. (4.51) 106 trabajo fin de máster La función objetivo 𝑧(𝑟)para este problema es complicada. Se conoce su valor sólo implícitamente como la solución de un problema de optimización en las variables de flujo 𝑥𝑘 𝑖𝑗. Además, para cualquier valor fijo de las variables de recursos 𝑟𝑘 𝑖𝑗, el problema de la directiva de recursos se descompone en un subproblema de flujo de red separado para cada servicio. Es decir, 𝑧(𝑟) = ∑𝑘∈𝐾𝑧𝑘(𝑟𝑘)con el valor 𝑧𝑘(𝑟𝑘)del 𝑘-ésimo subproblema dado por 𝑧𝑘(𝑟𝑘)=m ın ∑ 1≤𝑘≤𝐾 𝑐𝑘𝑥𝑘(4.52) s. a 𝑥𝑘=𝑏𝑘para todo 𝑘= 1,2,…, 𝐾. (4.53) 0≤𝑥𝑘 𝑖𝑗 ≤𝑟𝑘 𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴y para todo 𝑘= 1,2,…, 𝐾. (4.54) Propiedad 4.3.El problema de la directiva de recursos (4.45)-(4.48) es equivalente al problema de asignación de recursos (4.49)-(4.51) en el sentido de que 1. si (𝑥, 𝑟)es factible en el problema de la directiva de recursos, entonces 𝑟es factible en el problema de la asignación de recursos y 𝑧(𝑟)≤𝑐𝑥, y 2. si 𝑟es factible en el problema de asignación de recursos, entonces para algún vector 𝑥, (𝑥, 𝑟)es factible en el problema original y 𝑐𝑥 =𝑧(𝑟). Las propiedades (4.45) y (4.46) implican que en lugar de resolver el problema del flujo multiservicio directamente, se puede descomponer en un problema de asignación de recursos con una estructura de restricciones muy simple con una única restricción de desigualdad pero con una función objetivo complicada 𝑧(𝑟). Para encontrar 𝑧(𝑟) para cualquier elección del vector de asignación de recursos 𝑟, se necesita resolver 𝐾 problemas de flujo monoservicio. Otra forma de ver la función objetivo 𝑧(𝑟)es como el coste del problema lineal (4.45)- (4.48) en función de los parámetros 𝑟de la parte derecha. Es decir, cualquier valor 𝑟 para el vector de asignación define los valores de los parámetros del lado derecho para este problema lineal. Propiedad 4.4.Sea 𝑅el conjunto de asignaciones para el cual el problema lineal minimiza {𝑐𝑥 ∶𝑥=𝑏, 0≤𝑥≤𝑟}. Sea 𝑧(𝑟)el valor de este problema lineal en función del parámetro 𝑟del lado derecho. En el conjunto 𝑅, la función objetivo 𝑧(𝑟) es una función convexa, lineal a trozos de 𝑟. 4. problema de flujo multiservicio 107 4.8 Resolver los modelos de directiva de recursos Dado que la función 𝑧(𝑟)no es diferenciable (porque es lineal a trozos), no se puede usar métodos de gradiente para resolver el problema de asignación de recursos. En cambio, se podría utilizar otros métodos. Por ejemplo, podría buscarse una mejora local en 𝑧(𝑟)utilizando un método heurístico. Podría usarse un ”método de arco por vez” sumando 1 a 𝑟𝑘′ 𝑝𝑞 y restando 1 de 𝑟𝑘′′ 𝑝𝑞 para un arco (𝑝, 𝑞)en dos servicios 𝑘′y𝑘′′ , eligiendo el arco y los servicios en cada paso en base a algún criterio (por ejemplo, las opciones que den la mayor disminución en el valor de la función objetivo en cada paso). Este método es fácil de implementar pero no asegura la convergencia a una solución óptima. Puede visualizarse como cambiar la asignación de recursos en cada paso utilizando la fórmula: 𝑟←𝑟+𝜃𝛾, (4.55) con una longitud de paso de 𝜃= 1 y una dirección de movimiento 𝛾= (𝛾𝑘 𝑖𝑗)dada por 𝛾𝑘′ 𝑝𝑞 = 1,𝛾𝑘′′ 𝑝𝑞 = −1, y 𝛾𝑘 𝑝𝑞 = 0 para todas las demás combinaciones de arcoservicio. Sin embargo, tomando ideas prestadas de la optimización subgradiente, se podría usar un planteamiento de optimización eligiendo la dirección del movimiento 𝛾como un subgradiente correspondiente a la asignación de recursos 𝑟(ver Figura 4.3). Un método sería buscar una dirección de movimiento 𝛾, y una longitud de paso 𝜃 que simultáneamente mantengan la viabilidad y aseguren la convergencia hacia una solución óptima 𝑟del problema de asignación de recursos (4.49)-(4.51). Se adoptará un método de dos pasos. Primero, se determinará una dirección de subgradiente y una longitud de paso que asegure la convergencia, siempre que no se imponga restricciones a las asignaciones. Si moverse a 𝑟+𝜃𝛾 da una solución inviable, se transforma 𝑟en un punto 𝑟′que mantenga la factibilidad y asegure la convergencia. Para aplicar esta metodología, se debe poder responder dos preguntas básicas: 1. ¿Cómo se puede encontrar un subgradiante 𝛾de 𝑧(𝑟)en un punto 𝑟dado? 2. ¿Cómo se podría transformar una asignación de recursos no factible 𝑟+𝜃𝛾 de modo que sea factible para el problema de elección de recursos, es decir, satisfaga las restricciones ∑1≤𝑘≤𝐾𝑟𝑘 𝑖𝑗 ≤𝑢𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴, y 0≤𝑟𝑘 𝑖𝑗 ≤𝑢𝑘 𝑖𝑗 para todo (𝑖, 𝑗) ∈ 𝐴y todo 𝑘= 1,2,…, 𝐾? 108 trabajo fin de máster Asignación del recurso Subgradiante Figura 4.3: Subgradiante de la función z(r) Para responder a 1, recuérdese que, como se muestra en la Figura 4.3, un subgradiente 𝛾de 𝑧(𝑟)en el punto 𝑟=𝑟es cualquier vector que satisface la condición 𝑧(𝑟)≥𝑧(𝑟) + 𝛾(𝑟−𝑟)para todo 𝑟= (𝑟1, 𝑟2,…, 𝑟𝐾)con 𝑟𝑘∈𝑅𝑘.(4.56) En esta expresión, 𝑅𝑘es el conjunto de asignaciones de recursos para el servicio 𝑘 para el cual el subproblema (4.49)-(4.51) es factible. La siguiente propiedad muestra que para encontrar cualquiera de esos subgradiantes, se puede trabajar con cada una de las funciones constituyentes 𝑧𝑘(𝑟𝑘)de forma independiente. Propiedad 4.5.Sea 𝛾𝑘para 𝑘= 1,2,…𝐾un subgradiante de 𝑧𝑘(𝑟𝑘)en el punto 𝑟𝐾. Luego 𝛾= (𝛾1, 𝛾2,…, 𝛾𝑘)es un subgradiante de 𝑧(𝑟)en 𝑟= (𝑟1, 𝑟2,…, 𝑟𝐾). Este resultado muestra que para obtener un subgradiente de 𝑧(𝑟), se puede encontrar un subgradiente de cada función 𝑧𝑘(𝑟𝑘), lo cual requiere información sobre la sensibilidad de la solución de un problema de flujo en red (17.9) con cambios en las cotas superiores 𝑟𝑘 𝑖𝑗 de los flujos en arco. 4. problema de flujo multiservicio 109 Teorema 4.6 (Condiciones de Optimalidad de Holgura Complementaria). Una solución factible 𝑥∗es una solución óptima del problema de flujo de mínimo coste si y sólo si para algún conjunto de potenciales de nodo 𝜋, los costes reducidos y los flujos satisfacen las siguientes condiciones de optimalidad para cada arco (𝑖, 𝑗) ∈ 𝐴: Si 𝑐𝜋 𝑖𝑗 >0,luego 𝑥∗ 𝑖𝑗 = 0.(4.57) Si 0< 𝑥∗ 𝑖𝑗 < 𝑢𝑖𝑗,luego 𝑐𝜋 𝑖𝑗 = 0.(4.58) Si 𝑐𝜋 𝑖𝑗 <0,luego 𝑥∗ 𝑖𝑗 =𝑢𝑖𝑗.(4.59) Propiedad 4.6.Considérese el problema de flujo de red 𝑧(𝑞) = m ın{𝑐𝑥 ∶𝑥= 𝑏, 0≤𝑥≤𝑞}, con un vector paramétrico 𝑞con cotas superiores no negativas en los flujos en arco. Supóngase que 𝑥∗es una solución óptima a este problema cuando 𝑞=𝑞∗, y que 𝑥∗junto con el vector de coste reducido 𝑐𝜋satisfacen las condiciones de optimalidad del flujo de mínimo coste del Teorema (4.6). Se definen los vectores 𝜇= (𝜇𝑖𝑗)estableciendo 𝜇𝑖𝑗 ={0si 𝑥∗ 𝑖𝑗 < 𝑞∗ 𝑖𝑗 𝑐𝜋 𝑖𝑗 si 𝑥∗ 𝑖𝑗 =𝑞∗ 𝑖𝑗 (4.60) Entonces, para cualquier vector no negativo 𝑞′para el cual el problema es factible, 𝑧(𝑞′)≥𝑧(𝑞∗) + 𝜇(𝑞′−𝑞∗). Es decir, 𝜇es un subgradiente de la función objetivo paramétrica 𝑧(𝜇)en 𝜇=𝜇∗. Para aplicar la Propiedad 4.56, se resuelve cada subproblema (4.52)-(4.54). Después de resolver el subproblema del 𝑘-ésimo servicio con 𝑞=𝑟𝑘, se establece 𝛾𝑘igual al vector 𝜇. Luego se usa la propiedad 4.56 y se establece 𝛾= (𝛾1, 𝛾2,…, 𝛾𝐾). Estos resultados muestran cómo encontrar un subgradiente utilizando cualquier algoritmo de flujo de coste mínimo que produzca costes reducidos 𝑐𝜋así como un vector de flujo óptimo. Una vez que se haya obtenido un subgradiante 𝛾de 𝑧(𝑟), podría usarse la fórmula de subgradiante 𝑟←𝑟+𝜃𝛾 para encontrar la nueva asignación de recursos 𝑟. Si el vector de asignación de recursos resultante 𝑟es factible (satisface las restricciones ∑1≤𝑘≤𝐾𝑟𝑘 𝑖𝑗 ≤𝑢𝑖𝑗), se mueve a este punto. Sin embargo, si la nueva asignación de recursos 𝑟no es factible, es necesario modificarla para asegurar la factibilidad, se mueve en su lugar a algún otro punto 𝑟. Un método que asegura que el algoritmo converge con la elección apropiada de tamaños de paso (como se discutió en el Sección 2.2.6) es elegir 𝑟como el punto factible que está lo más cerca posible de 𝑟en el sentido de minimizar la cantidad ∑1≤𝑘≤𝐾∑(𝑖,𝑗)∈𝐴(𝑟𝑘 𝑖𝑗 −𝑟𝑘 𝑖𝑗)2. 110 trabajo fin de máster 4.9 Partición base El algoritmo puede utilizar la estructura de árbol de expansión de la base de programación lineal para acelerar los cálculos matriciales. El problema del flujo multiservicio es un problema lineal que combina dos conjuntos de restricciones: 1. restricciones de conservación de flujo independientes 𝑥𝑘=𝑏𝑘para cada servicio 𝑘= 1,2,…, 𝐾, y 2. restricciones de capacidad que vinculan los servicios a través de capacidades de arco compartidas. Debido a esta subestructura de red subyacente, es planteable si también podría usarse alguna forma de estructura de árbol de expansión para implementar el método simplex para el problema de flujo multiservicio. Considérese cualquier base de la formulación de programación lineal (4.1)-(4.4) del problema de flujo multiservicio. Dado que cada columna de corresponde al flujo en un arco de algún servicio, las columnas de la base definen un conjunto de subgrafos, uno para cada servicio. Sea 𝑘la submatriz de correspondiente al servicio 𝑘. Obsérvese que cada matriz 𝑘debe contener una base 𝑘de la matriz de incidencia de arco-nodo 𝑘. Los subgrafos correspondientes a cualquier base son árboles de expansión junto con algunos arcos adicionales. Además, dado que las restricciones del problema de flujo multiservicio tienen forma de bloques, también la tiene la base . La Figura 4.4 muestra la estructura de la base 𝑘. En lugar de trabajar directamente con la matriz base , mediante un cambio en las variables, se transformará esta matriz en otra matriz ′para poder realizar los pasos del método simplex de manera mucho más eficiente utilizando la estructura de red subyacente del problema. Supóngase que se establecen las restricciones de capacidad como restricciones de igualdad (4.5) agregando la variable de holgura no negativa 𝑠𝑖𝑗 a la restricción de capacidad para el arco (𝑖, 𝑗). La variable 𝑠𝑖𝑗 denota la capacidad de flujo no utilizada en el arco (𝑖, 𝑗). Considérese el sistema de balance de flujo y restricciones de capacidad como base del problema de flujo multiservicio, que se escribe como 𝑥𝐵=𝑏. En esta expresión, 𝑏= (𝑢, 𝑏1, 𝑏2,…, 𝑏𝐾)es el vector de los coeficientes del lado derecho de las ecuaciones y 𝑥𝐵es el conjunto de variables correspondientes a las variables básicas. Estas variables incluyen tanto las variables de flujo 𝑥𝑘 𝑖𝑗 como las variables de holgura 𝑠𝑖𝑗 para las restricciones de capacidad. 4. problema de flujo multiservicio 111 Restricciones de capacidad Servicio 1 Servicio 2 Servicio 3 Servicio K Servicio 4 Base (árbol de expansión) para cada servicio Figura 4.4: Partición de una base de programación lineal del problema de flujo multiservicio: base original. Para simplificar el sistema 𝑥𝐵=𝑏, supóngase que se hace un cambio en las variables dejando 𝑦𝐵=𝑥𝐵para una matriz no singular seleccionada apropiadamente. Si se sustituye 𝑦𝐵por 𝑥𝐵en el sistema 𝑥𝐵=𝑏, y ′=, el sistema se convierte en (𝑦𝐵) = ()𝑦𝐵=′𝑦𝐵=𝑏. Supóngase que se puede hacer este cambio en las variables de modo que la matriz ′tenga la estructura especial que se muestra en la Figura 4.5; en esta estructura especial, solo las variables 𝑦que aparecen en las restricciones de balance de flujo para el servicio 𝑘son las correspondientes a la base 𝑘. En el sistema transformado, sea la matriz de las restricciones de capacidad que no corresponden a ninguna de las columnas en las bases 𝑘para 𝑘= 1,2,…, 𝐾. Se hace referencias a como el trabajo base. En la representación pictórica de la base ′en la Figura 4.5, se ha reorganizado las columnas de la base de modo que aquellas en la base de trabajo aparezcan primero. Dos observaciones más: 112 trabajo fin de máster Restricciones de capacidad Servicio 1 Servicio 2 Servicio 3 Servicio 4 Servicio K Base (árbol de expansión) para cada servicio Figura 4.5: Partición de una base de programación lineal del problema de flujo multiservicio: base después de cambio de variables y reordenamiento de columnas. 1. Dado que un cambio en las variables, que corresponde a las operaciones de columna en la matriz , no cambia que esta matriz no es singular, la matriz ′ tampoco es singular. Este hecho implica que la base de trabajo no es singular. 2. Dado que las operaciones de columna en la matriz no cambian la solución del sistema 𝜋=𝑐𝐵, el vector de los multiplicadores simplex determinados por el sistema 𝜋′=𝑐𝐵′es el mismo que los determinados por el sistema original 𝜋=𝑐𝐵. (Cuando se realizan las operaciones de columna, se cambia tanto  como 𝑐𝐵;𝑐𝐵′denota el resultado de las operaciones de columna en los costes 𝑐𝐵). Estas dos observaciones y el cambio en las variables dan todos los ingredientes necesarios para implementar el método simplex de partición de bases para el problema de flujo multiservicio. Recuérdese que el método simplex requiere de dos tipos de cálculos matriciales: 1. Resolver un sistema de la forma 𝑥𝐵=𝑏o𝑥𝐵=𝑗para alguna columna 𝑗 4. problema de flujo multiservicio 113 de la matriz de coeficientes del modelo de programación lineal (se hace este cálculo en cada paso antes de introducir la columna 𝑗en la base). 2. Resolver un sistema de la forma 𝜋=𝑐𝐵para encontrar el vector 𝜋de multiplicadores simplex. En lugar de hacer estos cálculos de la matriz directamente, se usará la matriz ′. Los sistemas 𝜋=𝑐𝐵y𝜋′=𝑐𝐵′tienen las mismas soluciones. Supóngase que se divide el vector de multiplicadores simplex como 𝜋= (𝜋0, 𝜋1,…, 𝜋𝐾). El vector 𝜋0contiene los multiplicadores simplex para las restricciones de capacidad y el vector 𝜋𝑘contiene los multiplicadores simplex correspondientes a las restricciones del balance de flujo. 𝑥𝑘=𝑏𝑘para el servicio 𝑘. Sea 𝑐𝑊el subvector de 𝑐𝐵correspondiente a las columnas de 𝑊. Se resuelve el sistema 𝜋′=𝑐𝐵′de manera eficiente mediante el siguiente procedimiento de sustitución hacia adelante. Primero se resuelve el subsistema 𝜋′= 𝑐𝐵′resolviendo un sistema de ecuaciones lineales. Si se sustituye estos valores en el sistema 𝜋′=𝑐𝐵′, el sistema se descompone en subsistemas separados para cada servicio 𝑘con coeficientes modificados 𝑐𝐵′para 𝑐𝐵′. Además, los cálculos para cada subsistema son exactamente los necesarios para encontrar los multiplicadores simplex de un problema de flujo de red de un sólo servicio. Para resolver un sistema de la forma 𝑥𝐵=𝑏, primero se resuelve el sistema relacionado ′𝑦𝐵=𝑏usando sustitución hacia adelante. Es decir, si se divide el vector 𝑦𝐵en 𝑦𝐵= (𝑦0, 𝑦1, 𝑦2,…, 𝑦𝐾)correspondiente a las columnas 𝑊 , 1,2,…,𝐾de ′, primero se resuelve los subsistemas independientes 𝑘𝑦𝑘=𝑏𝑘. Luego se sustituye los valores de las variables 𝑦1, 𝑦2,…, 𝑦𝐾en el sistema ′𝑦=𝑏, dejando un sistema reducido de la forma 𝑦0=𝑢′(aquí 𝑢′denota las capacidades de arco modificadas que se obtiene después de sustituir los valores de las variables de flujo 𝑦1, 𝑦2,…, 𝑦𝐾 en las restricciones de capacidad transformadas). Después de resolver, en este sistema de ecuaciones lineales, se usa la transformación entre 𝑥𝐵e𝑦𝐵(es decir, 𝑥𝐵=𝑦𝐵) para encontrar los valores de las variables 𝑥𝐵. Supóngase que se aplica el procedimiento a un problema con 𝑝restricciones de capacidad (𝑝podría ser igual al número 𝑚de arcos si todos tienen restricciones de capacidad, o podría ser mucho menor que 𝑚). En lugar de resolver sistemas de la forma 𝑥=𝑏 en la matriz completa que tiene dimensión (𝑝+𝐾)×(𝑝+𝐾), se hacen 𝐾cálculos en árbol y se resuelve un sistema 𝑝×𝑝más pequeño 𝑦0=𝑢′. De manera similar, en lugar de resolver un sistema grande (𝑝+𝐾)×(𝑝+𝐾)de la forma 𝜋=𝑐𝐵se resuelve un solo sistema 𝑦0=𝑢′de dimensiones 𝑝×𝑝y se hace 𝐾cálculos en árbol. El resto de los pasos del método simplex (tasar las columnas para encontrar la variable 120 trabajo fin de máster dos nodos virtuales (asociados a los mismos nodos físicos de la arista original). Específicamente, cada arista está vinculada al nodo virtual identificado con una dirección. De esta manera, las aristas paralelas al eje 𝑋vinculan los nodos virtuales 𝑋, las aristas paralelas al eje 𝑌vinculan los nodos virtuales 𝑌y las aristas paralelas al eje 𝑍vincula los nodos virtuales 𝑍. 3. Los nodos virtuales asociados a un mismo nodo físico están vinculados a través de aristas virtuales. Sus longitudes son cero (ya que enlazan nodos con las mismas coordenadas), pero tendrán un coste positivo que representa el uso de un codo en la ruta. En la Figura 5.1 se muestra una ilustración de la explosión de nodos. En la imagen de la izquierda se muestra un nodo 𝑣de un cruce, que está unido a otros seis nodos. En la imagen de la derecha se muestra los tres nodos explotados, 𝑣𝑋,𝑣𝑌y𝑣𝑍, del nodo 𝑣. El nodo 𝑣𝑋(respectivamente, 𝑣𝑌,𝑣𝑍) está unido sólo a los dos nodos adyacentes mediante aristas paralelas al eje 𝑋(respectivamente, eje 𝑌, eje 𝑍) y con los otros dos nodos virtuales asociados al mismo nodo físico. Como puede verse, las tres aristas virtuales que unen 𝑣𝑋con 𝑣𝑌,𝑣𝑋con 𝑣𝑍y𝑣𝑌con 𝑣𝑍están dibujados con líneas discontinuas en esta imagen. En la Figura 5.2 se ilustra el uso de este grafo para modelar diferentes tipos de giros en una ruta. En la imagen de la izquierda, para pasar de 𝑎 a𝑏no se utilizan codos y se utiliza un único nodo explotado (𝑣𝑍). Por otro lado, en la imagen central (respectivamente, derecha), girar desde una arista paralela al eje 𝑍 hacia una arista paralela al eje 𝑌(respectivamente, Eje 𝑋) requiere atravesar la arista virtual (𝑣𝑌;𝑣𝑍) (respectivamente, (𝑣𝑍;𝑣𝑋)) agregando un coste por usar ese codo. El resultado es un grafo que representa el espacio factible donde se enrutan las tuberías. Figura 5.1: Ilustración de nodos explotados y aristas virtuales. 5. métodos basados en los flujos en redes para el problema del óptimo enrutamiento de tuberías en diseños navales 121 Figura 5.2: Ilustración del uso de nodos explotados y aristas virtuales para modelar diferentes tipos de giros en un camino del nodo 𝑎al nodo 𝑏. Se denota por  𝐺= (𝑉;𝐸)el grafo construido como se ha descrito anteriormente, donde 𝑉es el conjunto de nodos (virtuales) y 𝐸el conjunto de aristas. El conjun- to 𝐸incluye el conjunto de aristas virtuales, denotadas por 𝐸𝑣, y las réplicas de las aristas físicas iniciales. Dado que  𝐺está enmarcado en ℝ3, cada nodo 𝑣∈𝑉también se identifica con sus coordenadas 3D, (𝑥𝑣, 𝑦𝑣, 𝑧𝑣) ∈ ℝ3. Se denota por 𝑑𝑣𝑣′= 𝑑((𝑥𝑣, 𝑦𝑣, 𝑧𝑣),(𝑥𝑣′, 𝑦𝑣′, 𝑧𝑣′)) la distancia euclidia entre los nodos 𝑣, 𝑣′∈𝑉y por 𝑑𝑒𝑒′ la distancia euclidia entre las aristas 𝑒, 𝑒′∈𝐸(que representa tanto un punto en ℝ3 si 𝑒∈𝐸𝑣o un segmento en ℝ3si 𝑒∈𝐸⧵𝐸𝑣). Se supone que se tiene un conjunto finito de servicios, 𝐾, y cada uno de ellos requiere su propia tubería para determinar el diseño de su ruta. Cada servicio 𝑘∈𝐾se define mediante un par (𝑠𝑘, 𝑡𝑘) donde 𝑠𝑘∈𝑉es la fuente del servicio y 𝑡𝑘∈𝑉es el destino del servicio. Dependiendo del servicio 𝑘∈𝐾, se define una estructura de costes, 𝑐𝑘 𝑒≥0para cada arista del grafo 𝑒∈𝐸. Estos costes dependen de muchos factores como, por ejemplo, la longitud de las aristas, el coste de un codo, las preferencias del diseñador al enrutar las tuberías, etc. En la Sección 5.4.2 se describirá una estructura de costes detallada para el caso de estudio que trata este trabajo. 5.3 Un modelo basado en el flujo de red multiservicio (Faltan referencias y dibujos) En esta sección se describe el modelo de programación matemática propuesto para resolver el problema de enrutamiento de tuberías que impone requisitos técnicos de diseño naval. En una primera aproximación, se puede modelar el problema como un 122 trabajo fin de máster modelo de flujo de red multiservicio de mínimo coste (PFRMC, para abreviar). Se denota con 𝐺= (𝑉 , 𝐴)(donde 𝐴= {(𝑖, 𝑗) ∪ (𝑗, 𝑖) ∶ 𝑒= {𝑖, 𝑗} ∈ 𝐸}es el conjunto de arcos) es la versión dirigida del grafo  𝐺descrita en el Sección 5.2 previa. Similar a lo que ocurre en la versión no dirigida, se asume que el conjunto 𝐴incluye el conjunto de arcos virtuales, denotado por 𝐴𝑣, y las réplicas de los arcos físicos y se denota con 𝑑𝑎𝑎′la distancia euclidea entre los arcos 𝑎, 𝑎′∈𝐴(obsérvese que 𝑑𝑎𝑎′= 0 si 𝑎y𝑎′son los arcos correspondientes a las dos direcciones opuestas de la misma arista). Además, supóngase que, para cada servicio 𝑘∈𝐾, el coste de cada arco en 𝐴es 𝑐𝑘 𝑖𝑗 =𝑐𝑘 𝑗𝑖 =𝑐𝑘 𝑒. El objetivo del PFRMC es enrutar conjuntamente todos los servicios de 𝐾desde sus respectivas fuentes, 𝑠𝑘, a sus destinos, 𝑡𝑘, a través de la red 𝐺a un mínimo coste. Ejemplo 5.1.En la Figura 5.3, a la izquierda, se muestra un escenario que bien podría ser parte de un barco. El escenario mide tiene 3 servicios y un obstáculo de dimensiones que podría ser una viga. También en la Figura 1, a la derecha, se muestra la discretización del escenario en un mallado. Este mallado representa una versión simplificado del grafo 𝐺= (𝑉 , 𝐴)que realmente se usa en la resolución del problema de flujo multiservicio puesto que en el mallado los nodos no están explotados. Sea 𝑖un nodo genérico del mallado, entonces 𝑖𝑥, 𝑖𝑦, 𝑖𝑧son sus correspondientes nodos virtuales y sea 𝑛el número de nodos en el mallado. De manera que 𝑉= {𝑖𝑥, 𝑖𝑦, 𝑖𝑧}con 𝑖= {1,…, 𝑛}. A su vez, de acuerdo con la descripción de grafo dirigido 𝐺,𝐸𝑣= {(𝑣𝑖𝑥, 𝑣𝑖𝑦),(𝑣𝑖𝑥, 𝑣𝑖𝑧),(𝑣𝑖𝑦, 𝑣𝑖𝑧)} son las aristas virtuales resultantes de explotar los nodos y 𝐸⧵𝐸𝑣= {(𝑣𝑖𝑥, 𝑣𝑖+1𝑥),(𝑣𝑖𝑦, 𝑣𝑖+1𝑦),(𝑣𝑖𝑧, 𝑣𝑖+1𝑧)} las aristas reales que aparecen en el mallado. Para describir el modelo de programación matemática para el PFRMC, se usa el siguiente conjunto de variables de decisión binarias: 𝑥𝑘 𝑖𝑗 ={1si la ruta del servicio 𝑘utiliza el arco (𝑖, 𝑗), 0otro caso 𝑛 > 0.(5.1) Con este conjunto de variables, el coste total de enrutar los servicios a través de la red se puede expresar como: ∑ (𝑖,𝑗)∈𝐴∑ 𝑘∈𝐾 𝑐𝑘 𝑖𝑗𝑥𝑘 𝑖𝑗.(5.2) El siguiente conjunto de restricciones lineales asegura la representación correcta del 5. métodos basados en los flujos en redes para el problema del óptimo enrutamiento de tuberías en diseños navales 123 Figura 5.3: A la izquierda, escenario con tres fuentes en azul, tres sumideros en rojo y un obstáculo que se interpone entre ellos. A la derecha, mallado resultante de discretizar el espacio del escenario a la izquierda. problema: ∑ 𝑘∈𝐾∑ 𝑗∈𝑉 𝑥𝑘 𝑖𝑗 ≤1,∀𝑖∈𝑉 , (5.3) ∑ 𝑗∈𝑉 𝑥𝑘 𝑖𝑗 −∑ 𝑗∈𝑉 𝑥𝑘 𝑗𝑖 = 0,∀𝑘∈𝐾, 𝑖 ∈𝑉(𝑖≠𝑠𝑘, 𝑡𝑘),(5.4) ∑ 𝑗∈𝑉 𝑥𝑘 𝑠𝑘,𝑗 = 1,∀𝑘∈𝐾, (5.5) ∑ 𝑖∈𝑉 𝑥𝑘 𝑖,𝑡𝑘= 1,∀𝑘∈𝐾, (5.6) 𝑥𝑘 𝑖𝑗 ∈ {0,1} ∀(𝑖, 𝑗) ∈ 𝐴, 𝑘 ∈𝐾. (5.7) donde (5.3) establece que por un arco sólo pase un servicio, (5.4) impone que de un nodo sale tanto flujo como entra en él, (5.5) establece que el flujo que sale de una fuente debe ocupar una arista adyacente, (5.6) impone que una sólo por una arista 124 trabajo fin de máster adyacente a un sumidero debe pasar flujo y (5.7) limita las variables del modelo a valores binarios. La solución del PFRMC es un conjunto de rutas no superpuestas (una para cada servicio) que conectan las fuentes con los sumideros. Se sabe que este problema es 𝑁𝑃 - duro, incluso para el caso de dos servicios (ver [2]). Además, en el enrutamiento de tuberías hay dos requisitos técnicos que deben considerarse para construir diseños factibles: una mínima distancia entre las tuberías y una mínima distancia entre dos codos consecutivos de una tubería. El primer requisito permite separar adecuadamente los diferentes servicios para ser utilizados por las tuberías así como para mantener espacio suficiente para la manipulación de maquinaria en los camarotes de los barcos. El segundo, es una limitación física del trazado de la tubería, ya que no hay suficiente espacio para que una tubería gire dos veces en caso de que los codos estén colocados demasiado cerca. Las soluciones constan de segmentos en ℝ3sin dimensionalidad y, por ello, estas consideraciones no se tienen en cuenta en el PFRMC. Por lo tanto, para imponer esos requisitos técnicos, se deben agregar al problema las siguientes restricciones para modelar adecuadamente el problema de trazado de tuberías: 1. Distancia entre servicios: Dado una arista real (no virtual) usado por el servicio 𝑘, se requiere una mínima distancia de separación respecto a otros servicios. Sea 𝑅𝑘el radio del cilindro de la tubería 𝑘∈𝐾y sea Δ𝑘la mínima distancia de seguridad desde la superficie de la tubería 𝑘hasta cualquier otro elemento. La mínima distancia permitida entre los servicios 𝑘y𝑘′es entonces 𝑅𝑘𝑘′=𝑅𝑘+𝑅𝑘′+ m ax{Δ𝑘Δ𝑘′}. Se asegura el cumplimiento de este requisito imponiendo el siguiente conjunto de restricciones: ∑ 𝑘′∈𝐾∶𝑘′≠𝑘∑ 𝑎′∈𝐴⧵𝐴𝑣 𝑑𝑎𝑎′<𝑅𝑘𝑘′ 𝑥𝑘′ 𝑎′≤𝑀(1 − 𝑥𝑘 𝑎),∀𝑎∈𝐴⧵𝐴𝑣, 𝑘 ∈𝐾, (5.8) es decir, si 𝑎es un arco en la ruta del servicio 𝑘, entonces, no puede haber ningún arco en la ruta de cualquier otro servicio 𝑘′≠𝑘a una distancia menor que la mínima requerida. Aquí, 𝑀es una constante suficientemente grande. 2. Test de codos: Se requiere de un camino factible del servicio para verificar una distancia mínima , 𝐷𝑘≥0, entre codos consecutivos. Dado que los codos en el grafo no dirigido se identifican con aristas virtuales y luego, en el grafo dirigido con ar- 5. métodos basados en los flujos en redes para el problema del óptimo enrutamiento de tuberías en diseños navales 125 cos virtuales, este requisito se puede incorporar al modelo como sigue: 𝑥𝑘 𝑎+𝑥𝑘 𝑎′≤1,∀𝑎, 𝑎′∈𝐴𝑣∶𝑑𝑎𝑎′≤𝐷𝑘, 𝑎 ≠𝑎′, 𝑘 ∈𝐾, (5.9) es decir, si dos arcos virtuales (codos) diferentes del mismo servicio 𝑘∈𝐾están a una distancia menor o igual que 𝐷𝑘, entonces solo uno de ellos puede estar activo. Con las consideraciones anteriores, el problema de trazado de tuberías en un barco (PTTB) se puede formular como: m ın ∑ (𝑖,𝑗)∈𝐴∑ 𝑘∈𝐾 𝑐𝑘 𝑖𝑗𝑥𝑘 𝑖𝑗.(5.10) s. a ∑ 𝑘∈𝐾∑ 𝑗∈𝑉 𝑥𝑘 𝑖𝑗 ≤1,∀𝑖∈𝑉 , (5.11) ∑ 𝑗∈𝑉 𝑥𝑘 𝑖𝑗 −∑ 𝑗∈𝑉 𝑥𝑘 𝑗𝑖 = 0,∀𝑘∈𝐾, 𝑖 ∈𝑉(𝑖≠𝑠𝑘, 𝑡𝑘), (5.12) ∑ 𝑗∈𝑉 𝑥𝑘 𝑠𝑘,𝑗 = 1,∀𝑘∈𝐾, (5.13) ∑ 𝑖∈𝑉 𝑥𝑘 𝑖,𝑡𝑘= 1,∀𝑘∈𝐾, (5.14) ∑ 𝑘′∈𝐾∶𝑘′≠𝑘∑ 𝑎′∈𝐴⧵𝐴𝑣 𝑑𝑎𝑎′<𝑅𝑘𝑘′ 𝑥𝑘′ 𝑎′≤𝑀(1 − 𝑥𝑘 𝑎),∀𝑎∈𝐴⧵𝐴𝑣, 𝑘 ∈𝐾, (5.15) 𝑥𝑘 𝑎+𝑥𝑘 𝑎′≤1,∀𝑎, 𝑎′∈𝐴𝑣∶𝑑𝑎𝑎′≤𝐷𝑘, 𝑎 ≠𝑎′, 𝑘 ∈𝐾, (5.16) 𝑥𝑘 𝑖𝑗 ∈ {0,1},∀(𝑖, 𝑗) ∈ 𝐴, 𝑘 ∈𝐾. (5.17) Ejemplo 5.2.Continuando con el ejemplo 5.1 anterior, se resuelve este en la Figura 5.3. A la izquierda se resuelve el problema de flujo multiservicio sin introducir en el modelo las restricciones correspondientes al test de codos y al test de revestimiento. A la derecha aparece la solución del modelo una vez introducidas las restricciones correspondientes al test de codos y revestimiento. 126 trabajo fin de máster Figura 5.4: Solucion sin aplicar test de codos o test de revestimientos a la izquierda y aplicándolos a la derecha. El (PTTB) es un problema de flujo de red multiservicio de mínimo coste con restricciones adicionales que imponen el cumplimiento de los requisitos técnicos de las tuberías en la solución. Aunque el problema se modela como un problema de programación lineal de números enteros, en situaciones del mundo real tiene un gran número de variables ((|𝐴||𝐾|)) y un gran número de restricciones. En particular, las restricciones (5.15) y (5.16) imponen una carga considerable sobre el modelo ((|𝐴|2|𝐾|)). 5.3.1 Un método de solución exacta para el (PTTB) Las grandes dimensiones del modelo hacen inviable el uso de un solucionador de programación entera mixta incluso para instancias de pequeño tamaño. Este inconveniente se supera utilizando un método de ramificación y corte (sección 2.1) que incorpora las restricciones que se requieren, en lugar de considerar inicialmente todas. Esta estrategia es una forma eficiente de alcanzar la solución exacta del problema al 5. métodos basados en los flujos en redes para el problema del óptimo enrutamiento de tuberías en diseños navales 127 resolver una formulación incompleta con sólo algunas de las restricciones del modelo, mientras que las restricciones restantes, se incorporan en el modelo según sea necesario. Aunque en el peor de los casos, el procedimiento puede necesitar todas las restricciones que no se incluyeron inicialmente, en la práctica, solo se agrega una pequeña cantidad de ellas, reduciendo considerablemente el tamaño del problema. Específicamente, considerando al principio el siguiente problema maestro relajado: m ın ∑ (𝑖,𝑗)∈𝐴∑ 𝑘∈𝐾 𝑐𝑘 𝑖𝑗𝑥𝑘 𝑖𝑗.(5.18) s. a ∑ 𝑘∈𝐾∑ 𝑗∈𝑉 𝑥𝑘 𝑖𝑗 ≤1,∀𝑖∈𝑉 , (5.19) ∑ 𝑗∈𝑉 𝑥𝑘 𝑖𝑗 −∑ 𝑗∈𝑉 𝑥𝑘 𝑗𝑖 = 0,∀𝑘∈𝐾, 𝑖 ∈𝑉(𝑖≠𝑠𝑘, 𝑡𝑘),(5.20) ∑ 𝑗∈𝑉 𝑥𝑘 𝑠𝑘,𝑗 = 1,∀𝑘∈𝐾, (5.21) ∑ 𝑖∈𝑉 𝑥𝑘 𝑖,𝑡𝑘= 1,∀𝑘∈𝐾, (5.22) (5.23) Este problema no es más que un modelo de flujo de red múltiservicios que no considera las dos familias de restricciones técnicas que se requieren para una solución factible. Así, al resolver 5.18-, pueden obtenerse soluciones que verifican algunas de las restricciones (5.8) y/o (5.9). Para separar las restricciones violadas, se aplica un procedimiento enumerativo. En caso de que una de las restricciones no se verifique, se agrega al grupo de restricciones y el problema se resuelve nuevamente (eliminando la solución obtenida previamente). Más concretamente, para una solución dada del problema (5.18)-(5.3.1), 𝑥se verifica si viola (5.8) y/o (5.9) de la siguiente manera: 1. La violación de las restricciones (5.8) se comprueba examinando la solución, si una arista de un servicio se encuentra demasiado cerca de la arista de otro servicio, se añade una nueva restricción en el grupo de restricciones del problema principal. 2. Para verificar la violación de cualquiera de las restricciones (5.9) se examina la solución servicio a servicio. En caso de encontrar dos aristas virtuales demasiado cerca entre sí dentro de un mismo servicio, se añade una nueva restricción en el grupo de restricciones del problema principal. 128 trabajo fin de máster 5.3.2 Heurística El modelo de trazado de tuberías (PTTB) presentado en la Sección 5.3 requiere resolver un problema de programación de enteros, que en situaciones del mundo real, incluso el método de ramificación y corte no es capaz de cargar y resolver debido a la gran cantidad de variables y restricciones que tiene. En esta sección se proporciona un algoritmos basado en reducir la dimensión del modelo de programación entera. Una instancia del problema de trazado de tuberías es un ortoedro donde se han eliminado algunos ortoedros (obstáculos) más pequeños. Hay zonas de dicha región que rara vez son ocupadas por las rutas de las soluciones (por todos o algunos de los servicios). Se plantea un modelo que hace uso de esto para reducir su número de variables. Se propone una estrategia de recorte para las instancias. El procedimiento propuesto es iterativo y se describe en el pseudocódigo del algoritmo 7. El algoritmo comienza resolviendo el problema de la ruta más corta para cada servicio y determina la mejor solución para cada servicio independientemente de los servicios restantes. Luego, en lugar de resolver (PTTB) en todo el grafo, se resuelve el problema en la unión de tubos de sección cuadrada creada al agrandar las trayectorias obtenidas en una distancia desde la arista hasta el punto medio de la sección cuadrada inicial dada 𝛿𝑘(no declarando las variables que indican arcos fuera de las regiones cilíndricas). En caso de que el problema sea factible, proporciona una solución de (PTTB). De lo contrario, se aumenta 𝛿𝑘y se repite el procedimiento hasta que se obtiene la viabilidad. Aunque en el peor de los casos el algoritmo requiere resolver la instancia original del problema (el grafo completo), en la práctica, permite resolver el problema con un número considerablemente menor de variables. 5. métodos basados en los flujos en redes para el problema del óptimo enrutamiento de tuberías en diseños navales 129 distancia inicial para trazar los servicios alrededor de la solución del problema de su ruta más corta. 𝛿𝑘←𝑅𝑘+ Δ𝑘+𝛿para 𝑘∈𝐾 𝑖𝑡 ←1 𝑠𝑡𝑜𝑝 ←0 para 𝑘∈𝐾hacer Resuelva el problema del camino más corto para el servicio 𝑘en todo el grafo 𝐺∶ {𝑐𝑎𝑚𝑖𝑛𝑜𝑘}𝑘∈𝐾 fin mientras 𝑠𝑡𝑜𝑝 = 0 hacer Resuelve el (PTTB)(𝑖𝑡), donde todas las variables 𝑥𝑘 𝑖𝑗 para un servicio 𝑘 dado que están fuera del tubo alrededor de la ruta 𝑘y el radio 𝑘se fijan a cero. si (PTTB)(𝑖𝑡) es factible entonces 𝑠𝑡𝑜𝑝 = 1 en otro caso Aumentar 𝛿𝑘para todo 𝑘∈𝐾 fin 𝑖𝑡 =𝑖𝑡 + 1 fin Algoritmo 7: Reducción del número de variables. 5.4 Experiencia computacional (provisional, habría que cambiarlo) A continuación, se exponen los resultados de algunos experimentos computacionales realizados, con el fin de comparar empíricamente los métodos exacto y heurístico propuestos. 5.4.1 Instancias aleatorias Las instancias del grafo 𝐺= (𝑉 , 𝐸)se generan de la siguiente manera. Se genera un cubo de longitud de arista igual a 128 unidades que contenga un mallado de 17×17×17 nodos reales (𝑑= 17). Se generan cinco grupos diferentes de 15 servicios cada uno con un mismo radio de 4 unidades y estableciendo una distancia de seguridad entre servicios de 1 unidad. Para cada servicio, se ha generado aleatoriamente un origen y un destino en los nodos de la cuadrícula con 𝑑= 17 pertenecientes a los planos 𝑧= 0 136 trabajo fin de máster Figura 5.8: Representación gráfica de las tuberías resultantes al aplicar la solución de la Figura 5.7 para el escenario Ghenova con 10 servicios. A la derecha se muestra la solución sin mostrar los obstáculos. La figura de la izquierda muestra la misma solución integrada dentro del conjunto de obstáculos del escenario. 6 Conclusión En este trabajo se ha sintetizado en la medida de lo posible alguna de las principales herramientas disponibles para enfrentar los problemas de flujo multiservicio. También se han aplicado algunas de ellas para modelizar y resolver el trazado de tuberías y cables en barcos minimizando costes. Para ello se ha dispuesto de información de escenarios reales, esto es, partes de barcos; gracias a la colaboración de la empresa GHENOVA. Otras herramienas, sin embargo, se han estudiado en vistas a posteriores trabajos. En primera instancia, se ha partido tratando los algoritmos de ramificación y corte (Sección 2.1), pues es, herramienta ineludible para poder introducir los requisitos de diseño propios del diseño naval. Estos requisitos son, el test de codos (5.9) y la distancia entre servicios y obstáculos o test de revestimiento (5.8). Debido a las enormes dimensiones de los escenarios con los que irremediablemente se trabaja en este tipo de problemas, existen una cantidad colosal de estas restricciones y la forma de obtener una solución que cumpla estas restricciones es mediante la ramificación y corte. Seguidamente se ha considerado la relajación lagrangiana (Sección 2.2) y la generación de columnas (Sección 2.3). Se considera la relajación lagrangiana como una forma de relajar los problemas eliminando alguna de sus restricciones y aplicando un peaje en la función objetivo. Y la generación de columnas como una metodología para lidiar con problemas con un enorme número de variables, en la que se relajan sucesivos problemas y se compara su resultado con el de la relajación lagrangiana. Ambas dos herramientas son fundamentales en el estudio de los problemas de flujo multiservicio (Capítulo 4). Tras esto, se lleva a cabo un estudio del problema de flujo de mínimo coste mediante el algoritmo simplex en redes (capítulo 3). El problema de flujo de mínimo coste tiene siempre al menos una solución óptima de árbol de expansión y el algoritmo simplex en redes explota esta propiedad del problema recorriendo soluciones de ár- 138 trabajo fin de máster bol de expansión, sustituyendo en cada paso un arco del árbol por otro que no lo es. La comprensión de los problemas de flujo de mínimo coste es necesario para, en el siguiente capítulo, trabajar en el problema de flujo multiservicio. En conclusión de una parte más teórica de este trabajo, se estudia el problema de flujo multiservicio en el que varios servicios comparten una misma red y compiten por el espacio en ella para minimizar el coste total de los servicios (capítulo 4). Aquí hace uso de las herramientas anteriormente consideradas en el trabajo. Como aplicación práctica, en el siguiente capítulo (el 5) se modeliza el trazado de servicios con mínimo coste en el contexto del diseño naval. Se elabora un modelo que resuelve a optimalidad el problema y un heurístico necesario para escenarios grandes (todo escenario real). Se hace un estudio computacional de los resultados del modelo exacto y el algoritmo heurístico, así como se resuelve mediante el heurístico desarrollado un escenario real proporcionado por la empresa GHENOVA. En posteriores trabajos, se desea aplicar los resultados teóricos del Capítulo 4 en la resolución del problema descrito en el Capítulo 5. También se desea resolver otros problemas en diseño naval relacionados con este en colaboración con GHENOVA. Bibliografía [1] C. A. Floudas, P. M. P. E. Encyclopedia of Optimization. Springer, 2009. [2] Even, S., I. A., and Shamir, A. On the complexity of time table and multicommodity flow problems. 16th Annual Symposium on Foundations of Computer Science (1975), 184–193. [3] R. K. Ahuja, T. L. Magnanti, J. B. O. Network Flows: Theory, Algorithms and applications. Prentice hall, 1993. [4] Wolsey, L. A. Integer Programming. Wiley-Interscience, 1998.