Full text
ESCUELA DE DOCTORADO INTERNACIONAL DE LA UNIVERSIDAD DE SANTIAGO DE COMPOSTELA TESIS DE DOCTORADO CONTRIBUCIONES AL MODELADO Y CÁLCULO DE SOLUCIONES EN PROBLEMAS DE INVESTIGACIÓN OPERATIVA Jorge Rodríguez Veiga PROGRAMA DE DOCTORADO EN ESTADÍSTICA E INVESTIGACIÓN OPERATIVA SANTIAGO DE COMPOSTELA 2021
DECLARACIÓN DEL AUTOR DE LA TESIS “Contribuciones al modelado y cálculo de soluciones en problemas de investigación operativa” D. Jorge Rodríguez Veiga Presento mi tesis, siguiendo el procedimiento adecuado al Reglamento, y declaro que: 1. La tesis abarca los resultados de la elaboración de mi trabajo. 2. En su caso, en la tesis se hace referencia a las colaboraciones que tuvo este trabajo. 3. La tesis es la versión definitiva presentada para su defensa y coincide con la versión enviada en formato electrónico. 4. Confirmo que la tesis no incurre en ningún tipo de plagio de otros autores ni de trabajos presentados por mí para la obtención de otros títulos. En Santiago de Compostela, 14 de junio de 2021 Fdo. Jorge Rodríguez Veiga
AUTORIZACIÓN DE LA DIRECTORA DE LA TESIS “Contribuciones al modelado y cálculo de soluciones en problemas de investigación operativa” Dña. Balbina Virginia Casas Méndez INFORMA: Que la presente tesis, se corresponde con el trabajo realizado por D. Jorge Rodríguez Veiga, bajo mi dirección, y autorizo su presentación, considerando que reúne los requisitos exigidos en el Reglamento de Estudios del Doctorado de la USC, y que como directora de esta no incurre en las causas de abstención establecidas en la Ley 40/2015. De acuerdo con lo indicado en el Reglamento de Estudios de Doctorado, declara también que la presente tesis de doctorado es idónea para ser defendida en base a la modalidad de COMPENDIO DE PUBLICACIONES, en los que la participación del doctorando fue decisiva para su elaboración y las publicaciones se ajustan al Plan de Investigación. En Santiago de Compostela, 14 de junio de 2021 Fdo. Balbina Virginia Casas Méndez
Agradecimientos En primer lugar me gustaría agradecer a mi directora de tesis, Balbina Casas Méndez, todo el trabajo y apoyo que he recibido por su parte a lo largo de estos años. Gracias a ella comencé a tener curiosidad por la teoría de juegos y la optimización, y esta curiosidad ha conseguido mantenerla intacta durante todos estos años. De su mano realicé el Proyecto Final de Máster y posteriormente me animó para la realización de esta tesis. Quiero agradecerle la confianza que depositó al contratarme en el marco del proyecto LUMES, lo que supuso mi primer contrato con una empresa. Por todo esto y mucho más, gracias. Al igual que a Balbina, me gustaría agradecerle a Julio González Díaz todo el trabajo realizado conmigo. Tuve la fortuna de poder disfrutar muchos años trabajando bajo su tutela y todo lo aprendido creo que de un modo u otro está reflejado también en esta tesis. Gracias por todo el conocimiento, por implicarme en nuevos retos y por haber permitido que mi formación en mi etapa como doctorando se hubiese enriquecido llevándome a numerosos cursos y congresos. Aprovecho para dar mi profundo agradecimiento a todos los coautores de los trabajos que están recogidos en esta tesis, sin ellos esto no habría sido posible. Gracias Ángel González, Balbina Casas, David Rodríguez, Guido Novoa, Iván Gómez, José Luís Sáiz y María José Ginzo. Aprovecho para agradecer especialmente el trabajo de María José Ginzo por toda la ayuda y consejos que me ofreció durante todos estos años. También agradecer al Instituto Tecnológico de Matemática Industrial (ITMATI) la oportunidad de trabajar todos estos años junto a ellos en proyectos de transferencia, en concreto en los proyectos LUMES y las coi
Jorge Rodríguez Veiga laboraciones que surgieron con la empresa Repsol. Todo el conocimiento que he adquirido a lo largo de todos estos años es en gran parte gracias a ellos. Quiero agradecer a todos los IPs de los proyectos en los que me vi involucrado por confiar en mi trabajo y enseñarme tantas cosas, ya no solo del mundo de las matemáticas. En este sentido, gracias a: Balbina Casas, Wenceslao González, Beatriz Pateiro, Julio González, Alfredo Bermúdez y Francisco José Pena. Pero además, este agradecimiento es para todas esas personas con las que compartí horas de trabajo y de tiempo libre en ITMATI, porque todas esas personas conforman una gran familia para mí. Muchísimas gracias en especial a Adolfo Núñez (Fito), Andrea Vilar, Daniel Alves, Diego Rodríguez, Gabriel Álvarez, Irene Llana, Javier López, Joaquín Ossorio, Juan Bedoya, Manuel Cremades, Manuel Fontenla, Marcos Raydan, Oana Chis, Patricio Reyes (Pato) y Rabi Ouallam. Quiero dar las gracias también a mis compañeros de máster Alejandro Saavedra e Iria Roca, ya que esta andanza la comenzamos juntos y todas esas tardes y noches de risas y de agobios con una cerveza en la mano jamás las olvidaré. Por último, no puedo acabar estas palabras sin darle las gracias a mi familia. Agradecer a mis padres, hermano, abuelos y pareja por apoyarme en todo momento y confiar siempre en mí. En especial a mi madre, mi abuela y mi pareja; los tres pilares fundamentales de mi vida. Mi madre, porque desde pequeño me hizo entender las matemáticas como el juego más apasionante del mundo, algún día espero poder transmitir todo ese amor y conocimiento a mis hijos. A mi abuela porque siempre fue como una segunda madre, y aunque ya no esté aquí, sé que se sentiría muy orgullosa de mí. Y por último a mi compañera de vida, Noa, por toda su comprensión y apoyo incondicional recibido durante esta etapa. ii
Resumen La presente tesis abarca cuestiones tanto de teoría de juegos como de programación lineal entera mixta, conectadas por un hilo común, la investigación sobre modelado de problemas reales y el cálculo eficiente de soluciones. Primeramente, se presenta el cálculo de dos índices de poder, el índice con configuración y el índice de Banzhaf-Coleman generalizado, para juegos de mayoría ponderada con configuración de coaliciones. El trabajo novedoso consiste en el empleo de las funciones generatrices para la obtención de estos índices demostrando matemáticamente su idoneidad. Además, se realiza una extensión de los algoritmos a la clase más amplia de juegos de mayoría ponderada múltiple. Se presentan ejemplos de la vida real que muestran el alcance del modelo considerado y los algoritmos introducidos. El resto del trabajo se centra en la resolución de algunos problemas surgidos de la colaboración con la empresa Babcock España, líder española en servicios aéreos de emergencia. Los problemas reflejan los requerimientos relativos a la selección y organización óptima de recursos para la contención de incendios forestales. En la tesis se da solución a tres problemas concretos: la Selección y Asignación temporal de Recursos para la Contención de un incendio forestal (SARC), la Asignación de Aeronaves a Rutas de Vuelo (AARV) y la Asignación de Aeronaves a Puntos de Repostaje (AAPR). Además, debido a la complejidad del problema SARC, se realiza un estudio de la aplicabilidad de distintas técnicas de descomposición para mejorar la eficiencia en la resolución del mismo.
Jorge Rodríguez Veiga modelos abstractos de programación matemática, creados en 1759 por el economista y médico francés Francois Quesnay, para representar el flujo de mercancías a lo largo del proceso de producción y consumo; la investigación del matemático Charles Babbage [10], sobre el costo del transporte y la clasificación del correo (realizada en la Uniform Penny Post de Inglaterra en 1840); o el estudio de la gestión de inventarios, con el concepto de cantidad económica de pedido, desarrollado por el ingeniero de producción Ford W. Harris en 1913, entre otros. Pero es a comienzos del siglo XX cuando surge la actual investigación operativa (McCloskey [54], Rajgopal [66]). En los años precedentes a la Segunda Guerra Mundial, los ejércitos de Gran Bretaña, primero, y de los Estados Unidos después, ven la necesidad de hacer frente a la utilización eficaz de nuevas armas y nuevos desafíos relacionados con el despliegue de radares, la lucha antisubmarina y, en general, la dirección de operaciones militares complejas. Científicos en el Reino Unido, incluyendo al físico Patrick Blackett, el genetista Cecil Gordon, el zoólogo Solly Zuckerman, el biólogo Conrad Hal Waddington, el químico Owen Wansbrough-Jones, el estadístico Frank Yates, el matemático Jacob Bronowski y el físico y matemático Freeman Dyson, y en los Estados Unidos, con el físico y matemático George Dantzig, buscaron herramientas para tomar mejores decisiones en áreas como la logística o en problemas como el diseño de horarios de adiestramiento. Durante la Segunda Guerra Mundial, Blackett propulsó la fundación del campo de estudio conocido como investigación de operaciones. Durante la guerra, se mostró disconforme con las tácticas del bombardeo estratégico, usando la investigación de operaciones para demostrar que no tenía los efectos que los comandantes militares pensaban. Es interesante recoger la frase que Blackett incluyó en su informe Scientist for Operational Level (1941) destacando la importancia de la investigación de operaciones: “Se ha realizado bastante esfuerzo científico hasta ahora en la producción de nuevos dispositivos pero muy poco en el uso adecuado y eficaz de lo que hemos producido”. También otros matemáticos como John von Neumann comenzaron, en paralelo, a trabajar sobre una parte de la investigación de operaciones que hoy se conoce como teoría de juegos y que se ocupa del estudio 2
CAPÍTULO 1. INTRODUCCIÓN matemático de las situaciones conflictivas. Von Neumann propuso aplicar el lenguaje de la teoría de juegos y la teoría del equilibrio general para el estudio de la economía. Su primera contribución significativa fue el teorema minimax en 1928 [58] con información perfecta y dos jugadores,1 el cual extendió en 1944 en el libro Theory of Games and Economic Behavior [83], escrito junto con el economista Oskar Morgenstern. Una vez terminada la guerra, muchos de los científicos continuaron desarrollando los métodos empleados durante los años previos. Un hito fundamental en el desarrollo de la investigación operativa fue la introducción por George Dantzig, en el año 1947, del método símplex. El método símplex es en términos de aplicación generalizada, uno de los más exitosos de todos los tiempos, ya que la programación lineal domina el mundo de la industria. Además, es considerado por SIAM (Sociedad de Matemática Aplicada e Industrial) uno de los 10 algoritmos más relevantes del siglo XX. Es por ello que así como Blackett está considerado el padre de la investigación operativa, y von Neumann el de la teoría de juegos, Dantzig es considerado el padre de la programación lineal. Se ha de señalar que la presente tesis, que se enmarca dentro de la investigación operativa, abarca cuestiones tanto de teoría de juegos como de programación lineal. Ambas conectadas por un hilo común, la búsqueda eficiente de soluciones en problemas de la vida real. Hoy en día, el uso de modelos de investigación de operaciones es cada vez más frecuente como herramienta de apoyo en la toma de decisiones. Este mayor uso se explica, principalmente, por una serie de causas interconectadas. Entre ellas, cabe mencionar una apertura y un conocimiento cada vez más profundo de estas metodologías en un abanico de disciplinas cada vez más amplio, la aparición de problemas cada vez más complejos que se desea resolver, la mayor disponibilidad de software y el desarrollo de nuevos y mejores algoritmos de resolución. En la actualidad, las herramientas de la investigación de operaciones se aplican ampliamente a los problemas en los negocios, la industria y la sociedad en general. Ejemplos de ello es el uso de la investigación de operaciones en la industria petroquímica, las compañías aéreas, las 1El teorema establece que en ciertos juegos de suma cero existe una estrategia que permite a ambos jugadores minimizar su máxima pérdida. 3
Jorge Rodríguez Veiga finanzas, la logística y las administraciones públicas. La construcción de herramientas que hacen uso de la investigación operativa suele ser realizada por un equipo de especialistas junto con el cliente. Además, independientemente de la naturaleza del problema, cabe distinguir las siguientes etapas: 1. Definición del problema. 2. Formulación del modelo matemático. 3. Resolución del modelo. 4. Validación del modelo. 5. Interpretación y puesta en práctica de la solución. La definición del problema implica establecer una descripción clara y precisa de la situación a la que nos enfrentamos. Es en esta fase en la que el cliente tiene un papel crucial por ser la persona conocedora del problema. Por norma general, suelen ser necesarias varias reuniones para definir claramente las decisiones que se deberán de tener en cuenta, el objetivo que nos permite valorar las distintas decisiones, y las restricciones o limitaciones que el problema presenta. Este proceso de definición es crucial pues afectará de manera significativa a las conclusiones del estudio. La segunda etapa es la formulación de un modelo matemático. En ella se deberá de construir un modelo que represente el problema definido por el cliente. Este será representado por un conjunto de ecuaciones y expresiones matemáticas: las decisiones, que serán representadas por variables; una función objetivo, que permitirá medir la calidad de las mismas, y un conjunto de restricciones para que el modelo tenga en cuenta las limitaciones establecidas por el cliente. Este modelo representará una aproximación abstracta de la realidad. Una vez formulado el modelo matemático, podremos proceder a su resolución. Solucionar el modelo consiste en encontrar los valores de las variables, con el propósito de optimizar, si es posible, o cuando menos mejorar la eficiencia o la efectividad del sistema dentro del marco de referencia que fijan los objetivos y las restricciones del problema. La selección del método de resolución depende de las características del modelo. Muchos de los procedimientos de resolución tienen la característica de ser iterativos, es decir buscan la solución en base a la repetición de 4
CAPÍTULO 1. INTRODUCCIÓN la misma regla analítica hasta llegar a ella, si la hay, o cuando menos a una aproximación. Cuando se completa la primera versión del modelo, es inevitable que contenga errores. Para validar el modelo se deben realizar numerosos estudios de simulación. Con ellos se podrán identificar posibles errores para su corrección, así como comprobar que los resultados obtenidos son los deseados. La etapa de validación del modelo y de la resolución es crucial para determinar si el modelo propuesto hace lo que se supone que debe hacer. La etapa final se inicia con la puesta en práctica de la solución. En este punto es crucial traducir la solución encontrada a instrucciones y operaciones comprensibles para los individuos que intervienen en la toma de decisiones. En la Figura 1.1 se muestra cómo se deben emplear los métodos de la investigación operativa una vez se ha realizado el proceso de validación del modelo. En este caso diferenciamos la resolución de un problema matemático correctamente formulado frente a la intuición de los agentes involucrados en un proceso de toma de decisiones. Este esquema viene a reforzar la idea de que aunque en muchas ocasiones la intuición pueda ser una herramienta adecuada en la toma de decisiones, si se apoya en la interpretación de los resultados arrojados por un modelo de investigación operativa, esta puede fortalecerse tomando así decisiones más eficientes. Problema Modelo Resultados Decisión Formulación Resolución Interpretación Intuición Figura 1.1: Esquema de uso de la investigación operativa. Se puede concluir, que la investigación operativa se centra en encontrar solución a problemas de decisión en la vida real con el fin de maximizar o minimizar una o varias funciones objetivo y sujeto a una 5
Jorge Rodríguez Veiga serie de condicionantes o restricciones. En función del tipo de problemas se puede establecer una primera clasificación dentro de la investigación operativa que ya se fue perfilando a lo largo de la introducción histórica: la teoría de juegos y la optimización matemática. Cabe señalar que mientras los problemas de optimización involucran a un decisor, los problemas de teoría de juegos corresponden a más de un decisor. Para finalizar se describe a continuación la distribución de la tesis. En el presente capítulo se explica el contexto del trabajo realizado, estableciendo las bases históricas y también metodológicas de la disciplina en la que se enmarca la presente tesis, la investigación operativa. El Capítulo 2 establece las hipótesis y objetivos que se tratarán en la tesis, diferenciando los objetivos e hipótesis que se abordarán en los trabajos presentados. En el Capítulo 3 se realiza una revisión de las herramientas empleadas para la resolución de los problemas planteados. En este caso, se establece el marco teórico necesario para la comprensión del trabajo sobre el cálculo de índices de poder [72] y las herramientas empleadas para el cálculo de soluciones en los problemas de gestión de recursos en la contención de incendios forestales (trabajos [68], [70], [71], además del documento de trabajo [73]). El Capítulo 4 presenta una discusión general de los índices de poder, enmarcándolos dentro de la teoría de juegos, y de la gestión de recursos en la contención de incendios forestales, presentando el contexto en el que se producen este tipo de problemas y estableciendo la relación entre los mismos. En el Capítulo 5 se presentan los trabajos publicados y sometidos a revisión que conforman la tesis por compendio de artículos. Los trabajos incluidos en este capítulo son: [72] (Sección 5.1), [70] (Sección 5.2), [73] (Sección 5.3), [71] (Sección 5.4) y [68] (Sección 5.5). En el Capítulo 6 se presentan las conclusiones de la tesis, diferenciando el trabajo del cálculo de índices de poder, Sección 6.1, de los trabajos para una gestión eficiente de los recursos empleados en la contención de un incendio forestal, Sección 6.2. Además, se indica el trabajo futuro que se puede realizar siguiendo las líneas de trabajo establecidas. Por último, en el Apéndice A se presenta una descripción de las técnicas de descomposición empleadas en la Sección 5.3. La inclusión de este apéndice busca facilitar la comprensión del trabajo presentado en la sección previamente mencionada. 6
Capítulo 2 Hipótesis y objetivos En este capítulo, se explican las hipótesis y objetivos generales y específicos que se pretenden alcanzar con esta tesis, indicando además en qué publicación o publicaciones se abordan. 2.1. Hipótesis y objetivos generales En la actualidad existen diversos problemas que se pueden atacar desde el punto de vista de la investigación de operaciones. Algunos de los problemas ampliamente trabajados son los de planificación, asignación, localización, reparto de costes... Este trabajo se centra en la obtención de soluciones para dos tipos de problemas. El primero, Sección 2.1.1, consiste en la obtención de índices de poder en problemas de mayoría ponderada con configuración de coaliciones. El segundo, Sección 2.1.2, reside en la gestión de recursos mediante el modelado y resolución de problemas relacionados con la contención de incendios forestales. 7
Jorge Rodríguez Veiga 2.1.1. Índices de poder La teoría de juegos se clasifica en dos grandes áreas: juegos cooperativos y juegos no cooperativos. El primer problema abordado en esta tesis se enmarca dentro de los juegos cooperativos, donde los jugadores disponen de mecanismos que les permiten adoptar acuerdos vinculantes con otros jugadores. Los agentes involucrados en un juego se denominan jugadores. Dentro de los juegos cooperativos, cuando el pago que se puede garantizar cada coalición (o subconjunto de jugadores) como consecuencia de la cooperación entre sus miembros, se puede repartir de cualquier forma entre los miembros de dicha coalición, se habla de un juego cooperativo con utilidad transferible (abreviadamente, juego TU). Una de las principales metas en el estudio de los juegos TU es la definición de reglas de distribución (denominadas valores) de los pagos derivados de la cooperación de todos los jugadores entre cada uno de ellos. Dos de los valores más importantes son el de Shapley [76] y el de Banzhaf [11] (nombrado en muchas ocasiones como Banzhaf-Coleman en el contexto de índices de poder debido a la similitud respecto al trabajo realizado por Coleman [24]). Un caso particular de los juegos TU son los juegos simples, que cobran gran relevancia debido a su aplicación a las ciencias sociales y políticas. Un juego simple se puede definir a partir de un conjunto finito de jugadores y un conjunto de coaliciones denominadas ganadoras. Además, dentro de este tipo de juegos, están los juegos de mayoría ponderada, donde en vez de hablar de valor se usa el término índice de poder (debido a su uso en modelos de órganos de toma de decisiones en los que los acuerdos se toman por votación). El interés de estos juegos se suele centrar en conocer el poder o influencia, que tiene cada uno de los jugadores dentro del juego, una vez se conoce el resultado de la votación. En este contexto el valor de Shapley se renombra como índice de poder de Shapley-Shubik [77]. Además, parece lógico que en muchas ocasiones, los jugadores involucrados en el juego, tengan preferencias por unirse a unos jugadores frente a otros por diversos motivos, como pueden ser ideologías políticas o ubicaciones geográficas. De este modo, Owen generaliza los valores de Shapley-Shubik y Banzhaf-Coleman para juegos donde existe una parti8
CAPÍTULO 2. HIPÓTESIS Y OBJETIVOS ción del conjunto de jugadores en uniones a priori, nombrándolos como índice de Owen [63] y Banzhaf-Owen [60] respectivamente para juegos simples con estructura de coalición o uniones a priori. Por último, al extender los resultados a juegos simples con configuración de coaliciones (donde existe un recubrimiento del conjunto de jugadores tal que un jugador puede pertenecer a más de una coalición), entonces al valor de Owen y al índice de poder de Banzhaf-Owen se les conoce como índice con configuración [2] (también conocido como Owen-Shapley con configuración de coaliciones u Owen-Shapley CC) e índice de Banzhaf-Coleman generalizado [3] (también nombrado como Owen-Banzhaf con configuración de coaliciones u Owen-Banzhaf CC). Esta extensión a los juegos con configuración de coaliciones modela mejor determinadas situaciones reales en las que un jugador puede tener preferencias a colaborar con más de un grupo de jugadores. Por ejemplo, considérense las relaciones diplomáticas entre países. En la vida real, los países están organizados en coaliciones internacionales, no necesariamente disjuntas. Por ejemplo, Francia y España, entre otros, pertenecen a la Unión Europea y la Organización del Tratado del Atlántico Norte (OTAN), por otro lado Estados Unidos pertenece a la OTAN y al Tratado de Libre Comercio de América del Norte (TLCAN), mientras que México únicamente pertenece al TLCAN. Todos ellos pertenecen al Fondo Monetario Internacional (FMI). Dada la complejidad computacional para la obtención de estos índices de poder, en la actualidad un gran número de investigadores está trabajando en alternativas para su cálculo. Algunas de las más remarcadas son mediante el empleo de las funciones generatrices, técnicas de programación dinámica, métodos de enumeración o los métodos de Monte Carlo [53]. 2.1.2. Gestión de recursos en la contención de incendios forestales El segundo grupo de problemas pertenecen al campo de la optimización lineal entera mixta. En la actualidad, muchos de los problemas de logística, tanto en el ámbito público como privado, son abordados desde 9
Jorge Rodríguez Veiga el punto de vista de la optimización. Esto es debido a la versatilidad que tiene esta disciplina a la hora de añadir restricciones importantes en la toma de decisiones y al creciente deseo de minimizar costes, tiempos... La historia de la lucha contra los incendios forestales por parte de las administraciones tiene sus orígenes a finales del siglo XIX, con la aparición de los primeros parques nacionales en Estados Unidos. Tras establecerse el Parque Nacional de Yellowstone en 1872 como el primer parque nacional del mundo, la administración del parque asignó al ejercito de Estados Unidos la responsabilidad de su protección. Con su llegada al parque, detectaron que se producían gran cantidad de incendios en él y que se debía de decir cómo proceder a su control al no haber soldados suficientes. En ese momento se creó una política de extinción de incendios que posteriormente se aplicaría al resto de parques naturales. Varios eventos catastróficos en los Estados Unidos, como los incendios de Peshtigo, el del Cañón de Santiago y especialmente el Gran Incendio de 1910 contribuyeron a la creencia de que el incendio era un gran peligro contra el que había que luchar. La repetición de sucesos de gravedad, como los anteriormente mencionados a lo largo del siglo XX, convencieron a las autoridades de la necesidad de diseñar estrategias para controlar y extinguir los incendios, así como priorizar la seguridad y protección de los recursos de extinción (brigadas, vehículos contra incendios, helicópteros y aviones de extinción...). Es por ello que se incrementó la búsqueda de operaciones eficientes que garantizasen la seguridad de las vidas humanas (establecer zonas de seguridad, rutas de escape...). En la actualidad, debido a los nuevos avances tecnológicos que proporcionan ordenadores más potentes y redes de comunicación más rápidas y seguras, se ve la necesidad de realizar una gestión más eficiente de los recursos involucrados en la extinción de incendios forestales mediante técnicas de optimización y simulación. 10
CAPÍTULO 2. HIPÓTESIS Y OBJETIVOS 2.2. Hipótesis y objetivos específicos Los trabajos presentados en la tesis abordarán una serie de objetivos específicos. Estos perseguirán facilitar la toma de decisiones en problemas del mundo real mediante el planteamiento y resolución de los distintos problemas propuestos. A continuación se procederá a enunciar los objetivos específicos de la tesis y los trabajos en los que se acometen. El trabajo presentado en la Sección 5.1 tiene por objetivo el cálculo eficiente de los índices de poder de Banzhaf-Coleman y Owen para juegos de mayoría ponderada con configuración de coaliciones. Además se extiende este objetivo a una clase más amplia de juegos, los juegos de mayoría ponderada múltiple con configuración de coaliciones. Para ello, se presentan cuatro algoritmos que hacen uso de las funciones generatrices, realizando la programación de los mismos en el lenguaje abierto R e ilustrando el procedimiento propuesto con un ejemplo de la vida real tomado de las ciencias sociales. Los siguientes trabajos persiguen la búsqueda de soluciones en problemas de gestión de recursos en la contención de incendios forestales. Dentro de este campo, se presentan distintos problemas de gestión. El primero de ellos consiste en la Selección y Asignación temporal de Recursos para la Contención de un incendio forestal (SARC). Este problema tiene como objetivo facilitar la planificación del uso de los recursos involucrados en la contención de un incendio. Todo ello teniendo en cuenta el cumplimiento de la normativa española de no negligencia de frentes y periodos de descanso para pilotos y brigadas. Para abordar este objetivo se presenta en la Sección 5.2 un modelo de programación lineal entera para ayudar en la toma de decisiones. El modelo propuesto es programado utilizando el lenguaje R, facilitando su uso mediante la creación de una interfaz, la cual facilita la inclusión de datos, resolución y análisis de los resultados. Además, se realiza un estudio de simulación para analizar la viabilidad del modelo propuesto y se ilustra el uso del mismo mediante la definición de un ejemplo inspirado en una situación real. Pese a que el estudio de simulación realizado en la Sección 5.2 arroja buenos tiempos de resolución, se consideró importante mejorarlos. En la actualidad, el tiempo empleado para la toma de decisiones en el contexto expuesto es sumamente importante, pues una respuesta lenta podría de11
Jorge Rodríguez Veiga En los juegos TU, de forma general, parece natural exigir la propiedad de eficiencia, ya que es una forma adecuada de distribuir los costes o beneficios. Sin embargo, cuando se habla de juegos de votación, no se piensa en cómo se debe repartir el poder, sino que el interés radica en un índice que indique el poder que tiene cada jugador dentro del órgano de votación. Es por ello que las soluciones puntuales en los juegos simples y por ende de mayoría ponderada se denominan índices de poder. Definición 3.10. Un índice de poder sobre SI(N), con n=|N|, es una aplicación f:SI(N)−→ Rn que a cada juego (N, v)∈SI(N)le asigna un vector de Rncuya componente i-ésima representa el poder que se le asigna al jugador i. Aunque los juegos simples están dentro de la clase de los juegos TU, no todas las propiedades que se definieron sobre las soluciones tienen sentido dentro de esta nueva clase de juegos. Por ejemplo, la suma de dos juegos simples no da lugar a un juego simple, por lo que la aditividad no tiene cabida en el contexto de estos juegos. Es por ello que surge la propiedad de transferencia o aditividad para juegos simples. La propiedad de transferencia se apoya en la definición de dos juegos simples, la intersección de juegos simples (Definición 3.5) y la unión de juegos simples, la cual se define a continuación. Definición 3.11. Dados mjuegos simples (N, v1),...,(N, vm)∈SI(N), el juego simple (N, v1∨ · · · ∨ vm)viene dado para todo S⊂Npor, (v1∨ · · · ∨ vm)(S) := m´ax{v1(S),· · · , vm(S)}= =1si ∃j= 1, . . . , m :vj(S) = 1, 0en otro caso. Transferencia: Un índice de poder fsatisface la propiedad de transferencia si para todo par de juegos (N, v1),(N, v2)∈SI(N)se verifica que f(N, v1) + f(N, v2) = f(N, v1∨v2) + f(N, v1∧v2). Con la propiedad de transferencia, introducida por Dubey en 1975 [31], el propio Dubey en dicho artículo realizó una caracterización del índice de Shapley-Shubik y posteriormente, junto a Shapley, del índice de 18
CAPÍTULO 3. HERRAMIENTAS METODOLÓGICAS Banzhaf-Coleman [32]. Además, como en juegos simples v(S) = 1 si S∈W(N, v)(donde W(N, v)era el conjunto de coaliciones ganadoras, Definición 3.3) y v(S) = 0 en otro caso. Se tiene que v(S∪ {i})−v(S) = 1 si el jugador i hace ganador a S. Con esta idea surge el concepto de swing. Definición 3.12. Sea un juego (N, v)∈SI(N)y un jugador i∈N. Un swing para ies una coalición S⊂N\ {i}tal que S /∈W(N, v)y S∪ {i} ∈ W(N, v). Al conjunto de todos los swings del jugador i, se le denotará por S(i). Por tanto, con la idea de swing y con la propiedad de transferencia, se pueden reformular y caracterizar los índices de Shapley-Shubik y Banzhaf-Coleman de la siguiente forma. Teorema 3.13. El único índice de poder sobre SI(N)que verifica las propiedades de transferencia, jugador nulo, simetría y eficiencia es el índice de Shapley-Shubik. Dado un juego (N, v)∈SI(N)ei∈N, se define el índice de poder de Shapley-Shubik como, φi(N, v) := X S∈S(i) |S|!(n− |S| − 1)! n!. Teorema 3.14. El único índice de poder sobre SI(N)que verifica las propiedades de transferencia, jugador nulo, simetría y poder total es el índice de Banzhaf-Coleman. Dado un juego (N, v)∈SI(N)ei∈N, se define el índice de poder de Banzhaf-Coleman como, βi(N, v) := X S∈S(i) 1 2n−1. Juegos con uniones a priori y configuración de coaliciones En los modelos considerados hasta el momento se permitía realizar cualquier tipo de unión entre jugadores para formar coaliciones. Sin embargo, en la vida real esto no sucede, ya que debido a razones familiares, políticas o económicas unos jugadores pueden tener más afinidad para coaligarse con unos jugadores que con otros. Para modelar este tipo de 19
Jorge Rodríguez Veiga situaciones surgen los juegos TU con uniones a priori (también conocidos como juegos TU con estructura de coaliciones), donde la comunicación de los jugadores se ve restringida. Por ejemplo, en un parlamento los parlamentarios se agrupan de manera natural en partidos políticos distintos, los cuales a su vez se unen de acuerdo a sus ideologías. Definición 3.15. Un juego TU con un sistema de uniones es una terna (N, v, P)donde (N, v)∈TU(N)yP={P1, . . . , Pm}es una partición de N. Se denotará por U(N)al conjunto de juegos TU con uniones a priori y conjunto de jugadores N. En 1977, Owen propuso y caracterizó el valor de Shapley para juegos cooperativos con un sistema de uniones a priori. A este nuevo valor se le denomina valor de Owen [63]. Otras caracterizaciones del valor de Owen se pueden encontrar en [18], [84] o [20]. Además, Owen en 1981 también propuso la extensión del índice de Banzhaf-Coleman, que se denomina índice de Banzhaf-Owen [60]. Sin embargo, caracterizaciones de este nuevo valor no surgieron hasta los trabajos de Albizuri (2001) [1] y de Amer, Carreras y Giménez (2002) [8]. Para facilitar la escritura del valor de Owen y Banzhaf-Owen, dado un juego con uniones a priori (N, v, P)∈U(N)y un jugador i∈N, se denotará por Pi∈Pa la coalición a la que pertenece el jugador i. Definición 3.16. Dado un juego (N, v, P)∈U(N), se define el valor de Owen para todo jugador i∈Ncomo, φi(N, v, P) = X R⊂P\{Pi}X T⊂Pi\{i} g(Pi, R, T) (v(QRT ∪ {i})−v(QRT )) , siendo, QRT =[ Pk∈R Pk[T, y g(Pi, R, T) = |R|!(|P|−|R| − 1)! |P|! |T|!(|Pi|−|T| − 1)! |Pi|!. Definición 3.17. Dado un juego (N, v, P)∈U(N), se define el valor de 20
CAPÍTULO 3. HERRAMIENTAS METODOLÓGICAS Banzhaf-Owen para todo jugador i∈Ncomo, βi(N, v, P) = X R⊂P\{Pi}X T⊂Pi\{i} 1 2|P|+|Pi|−2(v(QRT ∪ {i})−v(QRT )) . Sin embargo los juegos con uniones a priori podrían no ser suficientes para representar situaciones de la vida real. En un gran número de ocasiones, un jugador no tiene por qué tener afinidad con un único grupo, sino que puede tener afinidad con más de uno por diversos motivos, tales como ideología política, preferencias geográficas... Es por ello que la clase de juegos con uniones a priori se puede ampliar a la clase de juegos con configuración de coaliciones, donde el conjunto de coaliciones no tiene por qué ser una partición del conjunto de jugadores sino un recubrimiento finito del mismo. Definición 3.18. Una configuración de coaliciones de Nes una colección finita de coaliciones de Ncuya unión es N, i.e., C={C1, . . . , Cm} es una configuración de coaliciones si, [ Ck∈C Ck=N. Definición 3.19. Un juego TU con configuración de coaliciones es una terna (N, v, C)donde (N, v)∈TU(N)yC={C1, . . . , Cm}es una configuración de coaliciones. Se denotará por CC(N)al conjunto de juegos TU con configuración de coaliciones y conjunto de jugadores N. En este contexto, se centrará el estudio en resultados para los juegos simples con configuración de coaliciones, denotando al conjunto de dichos juegos por SC(N). En 2006, Albizuri y Aurrekoetxea [3] y Albizuri, Aurrekoetxea y Zarzuelo [2] proponen el índice de Banzhaf-Coleman generalizado y el valor de configuración, para los juegos simples y juegos TU, respectivamente, con una configuración de coalición. Estas reglas generalizan el índice de poder de Banzhaf-Coleman y el valor de Owen, respectivamente. De forma similar a la presentada para el caso de juegos con uniones a priori, para facilitar la escritura del índice con configuración y del índice de Banzhaf-Coleman generalizado, dado un juego con una configuración 21
Jorge Rodríguez Veiga de coalición (N, v, C)∈CC(N)y un jugador i∈N, se denotará por Cial conjunto de coaliciones a los que pertenece el jugador i, i.e. Ci= {Ck∈C:i∈Ck}. Definición 3.20. Dado un juego (N, v, C)∈CC(N), se define el valor con configuración para todo jugador i∈Ncomo, φi(N, v, C) = X R⊂C\CiX Ck∈CiX T⊂Ck\{i} g(Ck, R, T)v(ˆ QRT ∪ {i})−v(ˆ QRT ), siendo, ˆ QRT =[ Cl∈R Cl[T, y g(Ck, R, T) = |R|!(|C|−|R| − 1)! |C|! |T|!(|Ck|−|T| − 1)! |Ck|!. Definición 3.21. Dado un juego (N, v, C)∈SC(N), se define el índice de Banzhaf-Coleman generalizado para todo jugador i∈Ncomo, βi(N, v, C) = X R⊂C\CiX Ck∈CiX T⊂Ck\{i} 1 2|C|+|Ck|−2v(ˆ QRT ∪ {i})−v(ˆ QRT ). En la Figura 3.1 se presenta un esquema de los resultados vistos hasta ahora en esta sección. En ella, se pueden distinguir tres columnas, en las que se posicionan los resultados pertenecientes a los juegos TU sin restricción en la comunicación, juegos con uniones a priori y juegos con configuración de coaliciones, respectivamente. En color gris se muestran los resultados relativos a los juegos TU, mientras que en color blanco se presentan los resultados en juegos simples. Por último, indicar que la forma de las cajas con esquinas redondeadas indica dónde se definió la solución, mientras que a las que les falta la esquina superior derecha muestran dónde se presentaron caracterizaciones de la solución asociada. En la Figura 3.1 se destacan los siguientes trabajos: Shapley (1953) [76], Owen (1977) [63], Albizuri, Aurrecoechea y Zarzuelo (2006) [2], Shapley y Shubik (1954) [77], Dubey (1975) [31], Banzhaf III (1965) [11], Coleman (1971) [24], Dubey y Shapley (1979) [32], Owen (1981) [60], Albizuri 22
CAPÍTULO 3. HERRAMIENTAS METODOLÓGICAS Shapley (Shapley, 1953) Banzhaf-Coleman (Banzhaf 1965 y Coleman 1971) Shapley-Shubik (1954) Banzhaf-Owen (Owen, 1981) Owen (Owen, 1977) Valor con configuración (Albizuri et al. 2006) Banzhaf-Coleman generalizado (Albizuri y Aurrekoetxea, 2006) Juegos sin restricción en la comunicación Uniones a priori Config. Coalic. Banzhaf (Owen, 1975) JUEGOS SIMPLESJUEGOS TU Dubey (1975) Shapley (1953) Dubey and Shapley (1979) Feltkamp (1995) Albizuri (2001), Amer et al. (2002) Owen (1977) Figura 3.1: Evolución histórica de los valores e índices de poder. (2001) [1], Amer, Carreras y Giménez (2002) [8], Albizuri y Aurrekoetxea (2006) [3], Owen (1975) [61] y Feltkamp (1995) [33]. El trabajo expuesto en la Sección 5.1 muestra resultados teóricos para el cálculo del índice con configuración y del índice de Banzhaf-Coleman generalizado en juegos de mayoría ponderada y mayoría ponderada múltiple con configuración de coaliciones. El cálculo de índices de poder en juegos simples de mayoría ponderada se puede simplificar mediante el uso de las funciones generatrices. Un método de análisis combinatorio que facilita el cálculo para contabilizar el número de elementos de un conjunto finito. Está técnica, pese a no poder emplearse en todo tipo de juegos, es ampliamente empleada en juegos simples de mayoría ponderada. Ejemplos de ello son el trabajo pionero de David Cantor [48] para el índice de poder de Shapley-Shubik, o los trabajos de Brams y Affuso [19] y Alonso-Meijide [6]. 23
Jorge Rodríguez Veiga Definición 3.22. Dada una sucesión {aj}j∈Nde números reales, la serie fa(x) = X j∈N ajxj, se denomina función generatriz de la sucesión {aj}j∈N, pudiendo ser finita o infinita. En la serie fa(x)la variable xno tiene significado propio y únicamente sirve para identificar cada ajcomo el coeficiente de xjen su desarrollo. Además, su uso se puede extender al empleo de distintos grupos de características, considerando en este caso, tantas variables como grupos de características haya. Si se suponen tres grupos de características, se puede plantear la función generatriz asociada a la sucesión {ajkl}j,k,l∈N como, fa(x, y, z) = X j∈NX k∈NX l∈N ajklxjykzl. Para analizar en detalle los fundamentos teóricos de las funciones generatrices se puede consultar el trabajo de Ríbnikov y Medkov [67]. Funciones generatrices en juegos de mayoría ponderada En la Sección 3.1 se enunció la definición de un juego de mayoría ponderada (Definición 3.4) y la definición de un swing (Definición 3.12). Para mejorar la notación, para un juego de mayoría ponderada (N, v)∈ SI(N)dado por [q;w1, . . . , wn], y para cualquier S⊂N, se denotará, w(S) := X i∈S wi. El índice de Shapley-Shubik para un juego simple de mayoría ponderada, (N, v)dado por [q;w1, . . . , wn], se puede reescribir en función del número de swings para el jugador i∈Nen coaliciones de tamaño r. Siguiendo la fórmula para el cálculo del índice de Shapley-Shubik para juegos simples presentada en el Teorema 3.13 se tiene que, φi(N, v) = X S∈S(i) |S|!(n− |S| − 1)! n!= n−1 X r=0 r!(n−r−1)! n!σi r(N, v), 24
CAPÍTULO 3. HERRAMIENTAS METODOLÓGICAS donde σi r(N, v)representa el número de swings para el jugador ique tienen cardinal ren el juego (N, v). El índice así expresado, se podrá calcular mediante el resultado expuesto en la siguiente proposición. Proposición 3.23. Sea (N, v)un juego simple de mayoría ponderada dado por [q;w1, . . . , wn]. El número de swings del jugador i∈Ncon cardinal res igual a σi r(N, v) = q−1 X k=q−wi νi k,r, siendo νi k,r =|{S⊂N:i /∈S, w(S) = kand |S|=r}|. El siguiente resultado, de David G. Cantor (se puede ver en [48]), proporciona la función generatriz de la sucesión {νi k,r}k,r∈N. Proposición 3.24. Sea (N, v)un juego simple de mayoría ponderada dado por [q;w1, . . . , wn]. La función generatriz de los números {νi k,r}k,r∈N definidos en la Proposición 3.23 viene dada por, Gφ i(x, z) = n Y j=1,j6=i (1 + xwjz). De la proposición anterior se deduce que, Gφ i(x, z) = n Y j=1,j6=i (1 + xwjz) = (1 + xw1z)· · · (1 + xwi−1z)(1 + xwi+1 z)· · · (1 + xwnz) =X S⊂N\{i}Y j∈S xwjz =X S⊂N\{i} xw(S)z|S| = n−1 X r=0 w(N\{i}) X k=0 νi k,rxk zr. Por lo tanto, los valores de νi k,r se identificarán como los coeficientes 25
Jorge Rodríguez Veiga de xkzral desarrollar Gφ i(x, z). Solo son de interés, debido a la Proposición 3.23, los coeficientes νi kasociados a xkcon q−wi≤k≤q−1. A continuación, se ilustrará el cálculo del índice de Shapley-Shubik haciendo uso de las funciones generatrices con un ejemplo. Ejemplo 3.25. Considérese N={1,2,3,4}como el conjunto de jugadores y un juego de mayoría ponderada representado por [3; 2,2,1,1]. Si se considera el jugador 1, siguiendo la Proposición 3.24, Gφ 1(x, z) = 4 Y j=2 (1 + xwjz) = (1 + x2z)(1 + xz)(1 + xz) =x4z3+ 2x3z2+x2z2+x2z+ 2xz + 1. Por la Proposición 3.23, los sumandos de interés son aquellos que tienen el exponente de xentre q−wi= 1 yq−1 = 2. Los sumandos a considerar son x2z2,x2zy2xz. Por lo tanto, teniendo en cuenta el grado de z, σ1 0(N, v) = 2 X k=1 ν1 k,0= 0, σ1 1(N, v) = 2 X k=1 ν1 k,1= 1 + 2 = 3,←− (x2z, 2xz) σ1 2(N, v) = 2 X k=1 ν1 k,2= 1,←− (x2z2) σ1 3(N, v) = 2 X k=1 ν1 k,3= 0. El índice de Shapley-Shubik para el jugador 1se calcula como φ1(N, v) = n−1 X r=0 r!(n−r−1)! n!σ1 r(N, v) = 1!2! 4! 3 + 2!1! 4! =3 12 +1 12 =1 3. El índice para los otros jugadores se obtiene de forma análoga. 26
CAPÍTULO 3. HERRAMIENTAS METODOLÓGICAS Para el cálculo del índice de Banzhaf-Coleman para un juego simple de mayoría ponderada, (N, v)dado por [q;w1, . . . , wn], el procedimiento es similar al seguido con el índice de Shapley-Shubik, salvo que en este caso no es necesario tener en cuenta el tamaño de las coaliciones. La siguiente proposición determina para un juego simple de mayoría ponderada el número de swings de un jugador. Proposición 3.26. Sea (N, v)un juego simple de mayoría ponderada dado por [q;w1, . . . , wn]. El número de swings del jugador i∈Nes igual a σi(N, v) = q−1 X k=q−wi νi k, siendo νi k=|{S⊂N:i /∈Sand w(S) = k}|. Teniendo en cuenta la proposición anterior, se puede reescribir el índice de Banzhaf-Coleman en función del número de swings (empleando la expresión presentada en el Teorema 3.14) del jugador i∈Ncomo, βi(N, v) = X S∈S(i) 1 2n−1=σi(N, v) 2n−1. Con estos resultados, Brams y Affuso [19], proporcionan el cálculo del índice de poder de Banzhaf-Coleman, para todo jugador, en un juego de mayoría ponderada mediante el uso de funciones generatrices. Esto lo realizan mediante la definición de la función generatriz que permite el cálculo de los números {νi k}k∈N. Proposición 3.27. Sea (N, v)un juego simple de mayoría ponderada dado por [q;w1, . . . , wn]. La función generatriz de los números {νi k}k∈N definidos en la Proposición 3.26 viene dada por, Gβ i(x) = n Y j=1,j6=i (1 + xwj). De la proposición anterior se puede deducir que, 27
Jorge Rodríguez Veiga Figura 4.1: Incendios con intervención de medios en 2018 en España. Copyright 2019 por Gobierno de España. Ministerio de Agricultura, Pesca y Alimentación [40]. Reimpreso con permiso. Pese a que el número de incendios y conatos ha ido decreciendo en la última década, la importancia de una buena gestión de los recursos es esencial. La magnitud del problema ocasiona el gasto de millones de euros por parte de las administraciones estatales, autonómicas y locales en la prevención y contención de los incendios forestales. Según los últimos informes del MAPA, el 44,4% de los siniestros en 2019 se han producido en el noroeste de la península, en la zona que abarca las comunidades de Galicia, Asturias, Cantabria y País Vasco y las provincias castellanas de León y Zamora. Además, aunque el número total de incendios al año se está reduciendo, existe una amenaza creciente, los Grandes Incendios Forestales (GIF).2En el 2018 pese a producirse tan solo 3 GIF, un 0,04 % sobre el 2Incendio en el que arden más de 500 hectáreas de superficie. 34
CAPÍTULO 4. DISCUSIÓN GENERAL total de siniestros ocurridos, estos supusieron el 20,97 % de la superficie afectada. Esta gran amenaza conlleva la necesidad de tomar decisiones de gran magnitud, teniendo en cuenta gran cantidad de variables que afectan a la toma de una decisión correcta. Por ejemplo, en Galicia, una de las comunidades autónomas con más incendios de España, entre el 13 y el 15 de octubre de 2017, se registraron más de 100 incendios forestales activos de forma simultánea, con un total de 45.000 hectáreas de terreno afectado por los incendios. El despliegue de medios realizado fue de 500 soldados, 35 brigadas, 220 motobombas y 20 recursos aéreos entre otros. Este ejemplo muestra la necesidad de hacer una buena gestión de los recursos disponibles debido a que aunque limitados, pueden ser numerosos, y su uso conlleva costes y riesgos que deben ser minimizados. El diseño de herramientas de apoyo a la gestión en logística es un campo ampliamente estudiado en investigación operativa. En la gestión de incendios forestales, destaca el estudio económico realizado por Headley [44] y Sparhawk [80] por ser precursores en este terreno. Desde un punto de vista más cercano a la investigación operativa, destaca el trabajo de Donovan y Rideout [29] por buscar una forma eficiente de gestionar los recursos disponibles en el ataque inicial a un incendio forestal. El objetivo del trabajo es conseguir la contención del incendio, minimizando el Cost Plus Net Value Change (C+NVC) [41], mediante la resolución de un problema de programación lineal entera mixta (MILP). El C+NVC considera los costes asociados al uso de los recursos y los costes ocasionados por las hectáreas de terreno afectadas por el incendio. Además, se menciona la importancia de considerar como costes por el efecto destructivo de los incendios, los asociados a la repoblación de los bosques, restauración de los bienes dañados... En la literatura relativa a la gestión en el contexto de incendios forestales, se distinguen tres problemas claramente diferenciados. El primero y segundo tratan sobre la prevención y la detección, respectivamente. El tercero y último, en el que se centra la presente tesis, se refiere a la gestión de los recursos que son necesarios para la contención de un gran incendio forestal [55]. 35
Jorge Rodríguez Veiga Ante este problema creciente, surge en el año 2010 en España el proyecto de investigación PROMETEO, con razón de mejorar la eficiencia en la lucha contra incendios. El proyecto liderado por Babcock España, hasta 2017 denominado Inaer, involucra tanto a entidades públicas como privadas y surge como una de las mayores apuestas en investigación aplicada en materia de lucha contra incendios. Contó con el trabajo de hasta 16 empresas y fue subvencionado en casi un 44% por el Centro para el Desarrollo Tecnológico Industrial (CDTI) mediante el programa de financiación para estimular la cooperación público-privada en investigación industrial, Programa CENIT, que contribuye a un mejor posicionamiento tecnológico del tejido industrial español. Como continuación del proyecto PROMETEO, surgen en los años 2013 y 2015 los proyectos LUMES y ENJAMBRE respectivamente. Ambos proyectos, promovidos por el Programa Estratégico de Consorcios de Investigación Empresarial Nacional (CIEN) y con la cofinanciación del CDTI, involucran a empresas tanto del sector público como privado. El objetivo de estos proyectos es desarrollar herramientas que faciliten y apoyen la rápida toma de decisiones seguras en la contención de incendios forestales. Uno de los propósitos de estos tres proyectos es mejorar la gestión de recursos en la contención de incendios forestales, pudiéndose clasificar las tareas de tal propósito de la siguiente manera: Diseño de un algoritmo para la estimación del perímetro por unidad de tiempo en un incendio forestal. Para esta tarea se hace uso de técnicas de estimación de conjuntos, a partir de imágenes térmicas obtenidas del incendio. Diseño de un algoritmo para la prevención de colisiones entre recursos aéreos en un incendio forestal. El objetivo es determinar un espacio seguro para la aeronave que permita prevenir situaciones de riesgo al intersecar dicho espacio con otros objetos. Para la realización de la tarea se ha de tener en cuenta la localización espacial de los recursos aéreos y la dirección y velocidad de desplazamiento. Con esta herramienta, se pretenden evitar situaciones de riesgo por colisión con otros recursos aéreos o con obstáculos difícilmente perceptibles por la vista, tales como tendidos eléctricos. Diseño de un algoritmo para el cálculo de la eficiencia de las descar36
CAPÍTULO 4. DISCUSIÓN GENERAL gas que realizan los recursos aéreos en la contención de incendios forestales. Para realizar la contención de un incendio forestal los recursos aéreos se abastecen de supresores o retardantes los cuales arrojan sobre los frentes del incendio. El objetivo de la tarea es estimar la eficiencia de las descargas, para posteriormente, determinar con mayor exactitud el número recursos necesarios para contener un incendio forestal. Diseño de un algoritmo para la gestión eficiente de recursos en la contención de un incendio forestal teniendo en cuenta la legislación vigente. El algoritmo ha de gestionar los recursos involucrados en la contención del incendio asegurando el cumplimiento de la normativa española de no negligencia de frentes (en cada frente ha de existir un número mínimo de recursos de contención) y periodos de descanso para pilotos y brigadas (establecidos en la Circular Operativa 16-B [79]). Aunque el autor de esta tesis trabajó en los distintos objetivos de estos proyectos, en el presente trabajo se detallará el último punto mencionado, relativo a una gestión eficiente de los recursos involucrados en la contención de un incendio forestal. En las siguientes secciones se diferencian tres problemas de gestión de gran relevancia en la contención de incendios forestales: la Selección y Asignación temporal de Recursos para la Contención de un incendio forestal (SARC), la Asignación de Aeronaves a Rutas de Vuelo (AARV) y la Asignación de Aeronaves a Puntos de Repostaje (AAPR). Selección y Asignación temporal de Recursos para la Contención de un incendio forestal (SARC) El problema SARC está inspirado en el problema propuesto por Donovan y Rideout [29]. Donovan y Rideout proponen un modelo MILP para la planificación de un ataque inicial a un incendio forestal con el objetivo de contener el incendio minimizando el C+NVC. La Figura 4.2 ilustra el tipo de solución que se busca resolviendo el problema de Donovan y Rideout. El modelo determina qué recursos han de ser seleccionados para trabajar en el incendio, y cuándo han de intervenir en su contención. En la figura se muestra que el incendio ha 37
Jorge Rodríguez Veiga sido detectado a las 9:00, y con los tramos representados en negro, se representa la actividad de los recursos en la contención del incendio. Además si la solución del problema es factible, se garantiza que cuando todos los recursos han abandonado el incendio es porque el incendio ha sido contenido. 09:00 10:00 Figura 4.2: Ilustración del problema de Donovan y Rideout. El problema SARC extiende el problema de Donovan y Rideout al considerar tres nuevos aspectos fundamentales. La primera necesidad consiste en incorporar políticas de uso de los recursos, teniendo en cuenta el tiempo máximo de trabajo diario, el tiempo máximo de vuelo sin la realización de descansos y el tiempo necesario para la consecución de un descanso. Estos requisitos surgen por la normativa vigente en España, 16-B, la cual establece jornadas de un máximo de 8 horas y 40 minutos de descanso cada 2 horas de vuelo para los recursos aéreos. La segunda, hace referencia a la necesidad de establecer un número mínimo y máximo de recursos de cada tipo. Esta restricción es necesaria para que los recursos trabajen en condiciones seguras, permitiendo establecer que ningún frente del incendio quede desatendido (empleando el número mínimo de recursos) y que, por ejemplo, el espacio aéreo del incendio no este saturado por aeronaves (número máximo de recursos aéreos). Por último, se adecuó el modelo del problema a situaciones distintas a la inicial. Se decidió extender el modelo para permitir la incorporación de estados iniciales de los recursos para poder resolver el modelo en cualquier situación. Este nuevo enfoque permite trabajar con la filosofía rolling horizon [75], que reacciona ante la evolución de la incertidumbre del problema, resolviendo un modelo actualizado cada vez que avanza el horizonte temporal de la optimización. 38
CAPÍTULO 4. DISCUSIÓN GENERAL La Figura 4.3 ilustra un ejemplo de solución para el problema SARC. El modelo determina, a partir de la situación actual del incendio, qué recursos han de ser seleccionados para trabajar en él, y cuándo han de intervenir en su contención siguiendo las restricciones anteriormente descritas. En la figura se ilustra cómo el modelo puede ser ejecutado en un periodo distinto al inicial (9:00), y cómo han de actuar los recursos seleccionados. Además de los tramos negros y continuos, ya descritos, se da la interpretación de otros dos. Los negros discontinuos representan los periodos de vuelo para desplazarse a bases, puntos de repostaje o al incendio. Por último, los tramos azules continuos representan periodos de descanso. 09:00 10:00 Figura 4.3: Ilustración del problema SARC. Para modelar el problema es necesario tener información relativa a los recursos disponibles y la evolución del incendio en una discretización del tiempo (periodos). Respecto a la información de los recursos, es necesario disponer de los costes de uso, la capacidad de contención, el tiempo necesario para llegar al incendio desde su ubicación, e información relativa a la situación inicial del recurso (en caso de estar realizando un descanso, cuánto tiempo ha descansado; cuánto tiempo ha trabajado sin realizar un descanso; cuánto tiempo de trabajo ha realizado en la jornada; así como información para saber si se encuentra sin trabajar, trabajando en el incendio actual, o trabajando en otro incendio) y a su política de uso (descansos, tiempos máximos de trabajo sin descansos...). Además, es necesario saber el número mínimo y máximo de cada tipo de recurso. Respecto a la información de la evolución del incendio, es necesario disponer de una estimación, para cada periodo de tiempo, del incremento del perímetro del incendio y el incremento de los costes asociados a la superficie de terreno afectada por el incendio. 39
Jorge Rodríguez Veiga La definición detallada del problema SARC se presenta en la Sección 5.2. Para su resolución se propone un modelo MILP que describe el problema por medio de la definición de conjuntos, parámetros y variables, así como por una función objetivo y un conjunto de restricciones. Además, en el artículo se ilustra el modelo con un ejemplo inspirado en una situación real, y se realiza un estudio de simulación para analizar los resultados en términos de función objetivo y tiempos de resolución, los cuales se consideran aceptables para horizontes temporales de 5 horas, considerando periodos de 10 minutos. Cabe destacar, que no se pudo hacer una comparativa entre las soluciones obtenidas por el modelo matemático propuesto y las decisiones tomadas en un caso real. Se solicitaron informes y datos históricos a distintas entidades, pero debido a problemas de confidencialidad, no se consiguió la información necesaria. Pese a no disponer de datos reales, para la construcción de los casos se recogió información de la prensa para establecer casos similares a los ocurridos en Galicia. Con esta información se establecieron situaciones realistas para la evolución del incendio, número de recursos que participan en la contención, características de los recursos... En la Sección 5.3, debido a la complejidad en la resolución del problema SARC, se realiza un estudio para mejorar los tiempos de resolución del problema. En este documento de trabajo se presentan los resultados obtenidos al resolver el problema mediante distintas técnicas de descomposición. Las técnicas empleadas son la descomposición de Benders [14], la descomposición branch and price [81] y la descomposición lagrangiana [43]. Además, se describe una forma alternativa para modelar el problema basándose en la técnica de descomposición de Benders junto con el conocimiento que se tienen del problema. En el documento de trabajo se realiza un estudio de simulación para estudiar la eficiencia de los métodos propuestos en términos de optimalidad y tiempos de resolución. En la Sección 5.5 se acompaña la definición del modelo matemático con la presentación de una interfaz gráfica diseñada mediante el paquete shinydashboard [21] de R. La interfaz facilita la inserción de datos en el modelo matemático, su resolución y la visualización de los resultados obtenidos para facilitar la interpretación de la solución del modelo. Además, se describe la necesidad de emplear modelos matemáticos para una 40
CAPÍTULO 4. DISCUSIÓN GENERAL actuación eficiente en la contención de un incendio forestal relacionando los tres problemas, el problema SARC y los dos que se describen a continuación (problemas AARV yAAPR). Asignación de Aeronaves a Rutas de Vuelo (AARV) En la revisión realizada por Martell (2015) sobre la gestión y toma de decisiones en incendios forestales [52] se indica que, “In the case of amphibious airtankers, the air attack officer must decide from which water body each airtanker will pick up water and when and where each airtanker will drop its load.” El problema de decisión descrito por Martell coincide con el problema AARV y será parte de la motivación, junto con las necesidades expuestas en los proyectos de investigación mencionados anteriormente, del trabajo expuesto en la Sección 5.4. En los incendios forestales los recursos aéreos se distribuyen en rutas de vuelo. Estas se definen como los circuitos elípticos que van desde un punto de abastecimiento de retardantes o supresores, hasta un frente del incendio, donde la carga de las aeronaves es arrojada. La distribución de los recursos aéreos en las rutas de vuelo es una tarea de gran relevancia para una contención eficiente del incendio. Esto es debido a que rutas de vuelo con gran distancia entre sus puntos de abastecimiento y descarga, ocasionará un menor número de descargas por hora de las aeronaves, empeorando por tanto el rendimiento de las mismas en las tareas de contención. Por tanto, para realizar una asignación eficiente de las aeronaves a las rutas de vuelo, se busca maximizar el rendimiento de la operación de extinción. Sin embargo, uno no se puede ceñir a seleccionar siempre las rutas de vuelo que proporcionen un mayor rendimiento. Esto es debido a que, en función del recorrido y tipología de los recursos aéreos, existe un número máximo de aeronaves por ruta. También es de gran importancia tener en cuenta el número de rutas de vuelo que comparten un mismo punto de abastecimiento; cómo se distribuye el porcentaje de retardantes o supresores en los distintos frentes del incendio (tarea 41
Jorge Rodríguez Veiga que deberá indicar el coordinador aéreo de medios aéreos en el caso de España), así como no desatender ninguno de sus frentes. En la Figura 4.4 se ilustra el problema AARV que se desea resolver. El modelo determina, la asignación de los recursos aéreos a las rutas de vuelo que ofrece un mejor rendimiento, asegurando el cumplimiento de las restricciones previamente indicadas. Por lo tanto se asegura la operación y se ajusta la importancia de los frentes al criterio de los coordinadores. Figura 4.4: Ilustración del problema AARV. Para modelar el problema, es necesario conocer los recursos aéreos seleccionados para la contención del incendio en un instante temporal, los puntos de abastecimiento de retardantes/supresores y los frentes del incendio que se desean atacar. Además, es necesario disponer de la capacidad de carga de los recursos aéreos, el número máximo de circuitos que pueden tener un punto de abastecimiento común, el porcentaje de retardantes/supresores que se desea en cada frente del incendio, el número máximo de recursos aéreos por ruta de vuelo y el número de descargas por hora en las rutas de vuelo. La definición detallada del problema AARV se presenta en la Sección 5.4. Para su resolución también se propone un modelo MILP y se ilustra el modelo propuesto con un ejemplo inspirado en una situación real. Además, se realiza un estudio de simulación para analizar los tiempos de computación que se obtienen en la resolución del problema ante distinto número de recursos aéreos, puntos de abastecimiento y frentes del incendio. Cabe señalar que las velocidades de crucero y las características to42
CAPÍTULO 4. DISCUSIÓN GENERAL pográficas no se consideran explícitamente. Sin embargo, el número de descargas por hora en las rutas de vuelo refleja implícitamente estos elementos. Concretamente, el número de descargas se estima mediante una recta de regresión, teniendo en cuenta las características de los recursos aéreos y las distancias entre los puntos de agua y los frentes. Esto permite al coordinador estimar el tiempo necesario para que una aeronave recorra la ruta, lo que permite calcular, el tiempo entre dos descargas de agua consecutivas. También es importante destacar, que aunque este problema es habitual en la contención de incendios forestales con recursos aéreos, no se encontró literatura al respecto para la realización de una comparativa con nuestro modelo. Asignación de Aeronaves a Puntos de Repostaje (AAPR) El combustible utilizado por los recursos aéreos es limitado, y en grandes incendios forestales, es necesario más de un punto de reabastecimiento de combustible para que las aeronaves puedan repostar. En España, tal y como se define en la Circular Operativa 16-B, el repostaje se realiza mientras los recursos aéreos realizan sus descansos en tierra. De acuerdo con estas regulaciones, las aeronaves deben tener un descanso mínimo de 40 minutos cada 2 horas de vuelo. Por lo tanto, para realizar operaciones eficientes, es importante que el repostaje no exceda el período de descanso establecido. Con la formulación del problema, se busca gestionar la asignación de recursos aéreos a las bases de repostaje, con el objetivo de minimizar el tiempo máximo de repostaje de las aeronaves. Es importante mencionar que una vez finalizada la tarea de repostaje, los recursos aéreos volverán a su plan de trabajo previamente asignado. Por lo tanto, es fundamental considerar el tiempo necesario para volar hasta el punto de repostaje, así como el tiempo necesario para regresar al incendio forestal. Además, se debe considerar la cantidad de recursos aéreos que pueden repostar simultáneamente en una base determinada. Por ejemplo, si un punto de repostaje es un camión cisterna con una sola manguera, el suministro simultáneo de combustible a múltiples recursos aéreos se hace imposible. De este modo, un recurso puede preferir esperar mientras otro completa 43
CAPÍTULO 5. TRABAJOS PUBLICADOS Y SOMETIDOS A REVISIÓN 5.1. Implementing generating functions to obtain power indices with coalition configuration Referencia del artículo J. Rodríguez-Veiga, G. I. Novoa-Flores y B. Casas-Méndez, «Implementing generating functions to obtain power indices with coalition configuration,» Discrete Applied Mathematics, vol. 214, págs. 1-15, 2016 Filiación autores: Jorge Rodríguez-Veiga1, Guido Ignacio Novoa-Flores2, Balbina Casas-Méndez3. Contribución: Implementación de los algoritmos y revisión de la redacción del artículo. Factor de impacto: 0.956. Categoría: Mathematics Applied-scie. Posición relativa: Nº118 de un total de 255 (Q2). Citas scopus: 4. Citas Google: 6. ISSN: 0166-218X. Enlace: https://www.sciencedirect.com 1Instituto Tecnológico de Matemática Industrial (ITMATI), España. 2Departamento de Estadística, Análisis Matemático y Optimización, Universidad de Santiago de Compostela, España. 3Grupo de Investigación MODESTYA, Departamento de Estadística, Análisis Matemático y Optimización, Universidad de Santiago de Compostela, España. 51
Jorge Rodríguez Veiga 5.2. An integer linear programming model to select and temporally allocate resources for fighting forest fires Referencia del artículo J. Rodríguez-Veiga, M. J. Ginzo-Villamayor y B. Casas-Méndez, «An integer linear programming model to select and temporally allocate resources for fighting forest fires,» Forests, vol. 9, n.o10, 583, págs. 1-18, 2018 Filiación autores: Jorge Rodríguez-Veiga1, María José GinzoVillamayor2, Balbina Casas-Méndez23. Contribución: Modelado, implementación, simulaciones y redacción del artículo. Factor de impacto: 2.116. Categoría: Forestry-scie. Posición relativa: Nº17 de un total de 67 (Q1). Citas scopus: 4. Citas Google: 10. ISSN: 1999-4907. Enlace: https://www.mdpi.com 1Instituto Tecnológico de Matemática Industrial (ITMATI), España. 2Grupo de Investigación MODESTYA, Departamento de Estadística, Análisis Matemático y Optimización, Universidad de Santiago de Compostela, España. 3Facultad de Matemáticas, Universidad de Santiago de Compostela, España. 52
CAPÍTULO 5. TRABAJOS PUBLICADOS Y SOMETIDOS A REVISIÓN 5.3. Application of decomposition techniques in a wildfire suppression optimization model Información del documento de trabajo J. Rodríguez-Veiga, D. Rodríguez-Penas, Á. M. González-Rueda y col., «Application of decomposition techniques in a wildfire suppresion optimization model,» inf. téc., 2021 Filiación autores: Jorge Rodríguez-Veiga12, David R. Penas2, Ángel M. González-Rueda3, María José Ginzo-Villamayor2. Contribución: Formulación e implementación de la descomposición de Benders así como de la reformulación del problema que consigue los mejores resultados. Además, diseño y ejecución de las simulaciones y redacción del artículo. 1Instituto Tecnológico de Matemática Industrial (ITMATI), España. 2Grupo de Investigación Modestya, Departamento de Estadística, Análisis Matemático y Optimización, Universidad de Santiago de Compostela, España. 3Grupo de Investigación Modes, Departamento de Matemáticas, Universidad de A Coruña, España. 53
Application of decomposition techniques in a wildfire suppression optimization model Jorge Rodríguez-Veigaa,c,∗, David R. Penasc, Ángel M. González-Ruedab, María José Ginzo-Villamayorc aTechnological Institute of Industrial Mathematics (ITMATI), Santiago de Compostela, Spain bModes Research Group, Department of Mathematics, University of A Coruña, A Coruña, Spain cModestya Research Group, Department of Statistics, Mathematical Analysis and Optimization, University of Santiago de Compostela, Santiago de Compostela, Spain Abstract Resource assignment modeling provides an automatic and fast decision support system for wildfire suppression logistics. However, this process generates challenging optimization problems in many real-world cases, and the computational time becomes a critical issue, especially in large-scale instances. Thus, to overcome that limitation, this work studies and applies a set of decomposition techniques such as augmented Lagrangian, branch and price, and Benders decompositions to a wildfire suppression model. Moreover, a reformulation strategy, inspired by Benders’ decomposition, is also introduced and demonstrated. Finally, a numerical study comparing the behavior of the proposals using different problem sizes is conducted. Keywords: Integer programming; assignment problems; wildfire management; decomposition techniques; Benders decomposition 1. Introduction Many critical problems in disaster management and logistics can be formulated and solved using the Operations Research (OR) framework (Van Wassenhove and Pedraza Martinez, 2012; Caunhye et al., 2012). Wildfires are a type of catastrophe with a high impact on the current world, in humanitarian, economic, and above all, ecological terms (European Commision, 2019). As their frequency and magnitude are growing at an alarming rate (Dennison et al., 2014; Úbeda ∗Corresponding author: [email protected] Preprint submitted to European Journal of Operational Research March 16, 2021
and Sarricolea, 2016), it is necessary to develop efficient methods to improve the prevention, detection, and planning in logistics related to wildfire suppression. Several research topics in the optimization of forest fire management have been proposed in recent years (Minas et al., 2012; Miller and Ager, 2013), resulting in new challenging problems from the OR point of view. A strategy for managing the resources involved in a wildfire is to use a rolling horizon methodology: a resource coordinator runs a mathematical programming model multiple times, updating the input data for each run according to the evolution of the fire (Rodríguez-Veiga et al., 2018a). For that reason, solving each optimization problem for each iteration cannot be an expensive task with the execution time being a critical issue. A possible solution to obtain planning as soon as possible is to design fast heuristic methods that obtain a feasible solution without optimal guarantees. However, when the objective is to reach the global optimum or at least a close enough solution, a good alternative is to apply mathematical decomposition methods (Conejo et al., 2006). Based on previous experience in resource management (Rodríguez-Veiga et al., 2018a,b), we applied a set of decomposition strategies to a Mixed-Integer Linear Programming (MILP) model focused on wildfire suppression. In detail, we propose some specific decomposition methods and a new model (based on an extension of one of the decompositions) to tackle this optimization problem with excellent results. The organization of this document is as follows. Section 1.1 covers the related work while Section 1.2 presents a brief revision of the optimization model to be addressed: a logistic scheduling model for wildfire suppression. Section 2 describes different ways to apply decomposition strategies to the previous model, and Section 3 proposes a customized method based on the Benders decomposition. The performance of our proposal is evaluated using simulated data of the wildfire model and solving instances of different sizes in Section 4. Finally, Section 5 summarizes the main conclusions of our study. 1.1. Related Work Several researchers have addressed logistic management in wildfire suppression through the OR methodology using well-known optimization models, such as vehicle routing problems (VRPs), facility location problems, or scheduling problems. In the case of the first, routing problems can help in different ways, either minimizing the distance to reach the fire point from the deposit or building evacuation paths. An example is proposed in Yang et al. (2019), where a multiobjective VRP model optimizes both the delivered cost and the travel time. Another issue could be to obtain the optimal location of forest fire attack facilities (operational bases where suppression resources are stored) while avoiding their overlap. These location problems are addressed in Lee (2006), where a 2
stochastic MILP approach proposes two different models: one model for the optimal location of these bases under uncertainty in the occurrence of fires and another model for planning the optimal deployment of firefighting. Scheduling problems in forest fire suppression are perhaps the most frequent topic in the literature and probably the most important topic in this context. In Martell (2007), the authors apply OR methodology to help decision-making issues in forest fire prevention, detection, deployment, and initial dispatch. Another example is Zhou and Erdogan (2019), where a stochastic model to manage resource assignation and resident evacuation is presented. As stated above, applying decomposition techniques simplifies optimization problems. In Conejo et al. (2006), several of the decomposition techniques are explained in detail, such as Benders’ decomposition, Dantzig-Wolfe’s decomposition, and Lagrangian decomposition. There are many examples in the literature about applying these techniques to MILP problems. We are especially interested in the problems that study logistics scheduling and contemplate recurring elements in the context of a wildfire, such as the existence of air resource management, regulations regarding the time limit for the use of resources, or the presence of temporal periods in the model. Benders’ decomposition algorithm is a popular strategy implemented with success in many works. For instance, Mercier et al. (2005) used it to combine aircraft routing and crew scheduling problems, seeking to obtain a solution where aircraft and crews are assigned to a flight. In Papadakos (2009), several airline scheduling optimization models were studied, and they were solved as an integrated model via an enhanced Benders decomposition method combined with accelerated column generation. Moreover, Romanski and Hentenryck (2016) consider prescriptive evacuation planning for a region threatened by a natural disaster, such as a flood, forest fire, or a hurricane. They propose a Benders decomposition that generalizes the two-step approach proposed in previous work for converging evacuation plans. Dantzig-Wolfe decomposition usually appears in the literature together with the branch and price algorithm (Barnhart et al., 1998) in the MILP context, obtaining excellent results in different case studies. Thus, works such as Martins et al. (2012) used these approaches to solve forestal harvest scheduling problems with constraints on the maximum clear-cut area. Another interesting example is Rios and Ross (2010). It applies a decomposition based on branch and price to an air-traffic scheduling problem, obtaining one flight per subproblem. Regarding the Lagrangian decomposition method, a good example of resource allocation appears in Nishi et al. (2005). They applied a decomposition approach to a problem related to the route planning of multiple automated guided vehicles, including avoiding collisions, with the goal of minimizing the transportation costs. As far as we know, there is no work dedicated to implementing mathematical 3
decomposition techniques in optimization models in the context of resource planning in wildfires. We hope to contribute to simplifying the application of these effective strategies in different forestry optimization problems. 1.2. Problem statement: a wildfire suppression model We begin our study from a wildfire extinguishment model, described in RodríguezVeiga et al. (2018a), where an MILP model is presented to find an optimal resource assignation to extinguish a wildfire. The resources are allocated in different periods to contain a wildfire, according to Spanish regulations for the nonnegligence of the fronts and rest periods for pilots and brigades. A short description of this model is shown in the following paragraphs. Table 1 shows a definition of the different sets and parameters used in the model, and Table 2 describes the decision and auxiliary variables. The objective function and the constraints of the model are formulated as follows: minimize X i∈I,t∈T Ci·uit +X i∈I Pi·zi+X t∈T NV Ct·yt−1+X g∈G,t∈T M0·µgt (1) subject to X t∈T P ERt·yt−1≤X i∈I,t∈T P Rit ·wit (2) ∀t∈ T , M ·yt≥X t0∈T t P ERt0·yt−1−X i∈I,t0∈T t P Rit0·wit0(3) ∀i∈ I, t ∈ T , Ai·wit ≤X t0∈T t trit0(4) ∀i∈ I :IT Wi= 1,si1+X t∈T2 (m+ 1) ·sit ≤m·zi(5) ∀i∈ I :IT Wi= 0,X t∈T sit ≤zi(6) ∀i∈ I, t ∈ T ,X t0∈T t t−T RPi+1 trit0≥T RPi·eit (7) ∀i∈ I, t ∈ T ,0≤crit ≤W Pi(8) ∀i∈ I, t ∈ T ,rit ≤X t0∈T t+RPi−1 t erit0(9) 4
Table 1: Sets and parameters of the problem. Sets Definition i,IIndex and set of resources. g,GIndex and set of different resources groups. t,TIndex and set of time periods with wildfire behaviour information {1, . . . , m}. i,IgIndex and set of resources associated with group g. t,Tt2 t1Index and set of time periods between the natural numbers t1 and t2,{max{1, t1}, . . . , min{m, t2}}. If the sub-index or superindex is omitted, it means that t1= 1 and t2=m,respectively. When t1is greater than t2, we will consider it an empty set. Parameters Definition CiCost per usage period of resource i∈ I. PiFixed cost for resource i∈ I selection. NV CtIncrease in the costs of the wildfire in period t∈ T . P ERtIncrement of the wildfire perimeter in period t∈ T . P Rit Performance of resource i∈ I in period t∈ T . IT WiValue 1 indicates resource i∈ I is currently working in this wildfire (0 otherwise). IOWiValue 1 indicates resource i∈ I is currently working in other wildfire (0 otherwise). CW PiNumber of periods used since the last final rest of i∈ I. CRPiNumber of rest periods used by resource i∈ I. CUPiNumber of usage periods in the day of resource i∈ I. AiNumber of periods needed by i∈ I to arrive to the wildfire. W PiMaximum allowed number of periods without breaks for i∈ I. T RPiNumber of periods needed by resource i∈ I to go from the resting point to the fire. RPiNumber of rest periods that resource i∈ I must do on a break before starting to work again. UPiMaximum number of allowed usage periods in a day for i∈ I. nMingt Minimum number of resources of group g∈ G working on the wildfire in the same period t∈ T . nMaxgt Maximum number of resources of group g∈ G working on the wildfire in the same period t∈ T . M0Positive constant that penalizes the breach of the minimum number of resources in each period. MPositive sufficiently large constant to establish wildfire containment. A suitable value is M:= Pt∈T P ERt. 5
Table 2: Variables of the problem. Variables Definition sit ∈ {0,1}It takes the value 1 if resource i∈ I is starting to be used in period t∈ T . trit ∈ {0,1}It takes the value 1 if resource i∈ I is travelling in period t∈ T . rit ∈ {0,1}It takes the value 1 if resource i∈ I is resting in period t∈ T . erit ∈ {0,1}It takes the value 1 if resource i∈ I is ending a break in period t∈ T. eit ∈ {0,1}It takes the value 1 if resource i∈ I is ending its work in period t∈ T . yt∈ {0,1}It takes the value 1 if the fire is not contained in period t∈ {0} ∪ T . µgt ∈NIt indicates the number of missing resources of group g∈ G to reach the corresponding minimum in period t∈ T . Aux Vars Definition uit Resource i∈ I has a task assigned in period t∈ T . uit := Pt0∈T tsit0−Pt0∈T t−1eit0. wit Resource i∈ I works fighting the wildfire in period t∈ T . wit := uit −rit −trit. ziResource i∈ I is selected to fight the wildfire. zi:= Pt∈T eit. criNumber of periods since the last ending break period of resource i∈ I in period t∈ T . crit := X t0∈T t (t+ 1 −t0)·sit0 −X t0∈T t (t−t0)·eit0−X t0∈T t rit0 −X t0∈T t W Pi·erit0 if IT Wi= 0 and IOWi= 0 (t+CW Pi−CRPi)·si1 +X t0∈T t 2 (t+ 1 −t0+W Pi)·sit0 −X t0∈T t (t−t0)·eit0−X t0∈T t rit0 −X t0∈T t W Pi·erit0 in the other case 6
weight given to each solution added to the master problem) to the BPD master problem. If the solution improves the objective function, then this new solution is added to the BPD master problem and solved again, giving new dual value solutions to update the BPD subproblem. Thus, this procedure is repeated until no more columns can be added. Note that the scheme we have just described must be applied at every node of a branch and bound procedure. 2.3. Benders decomposition Benders decomposition (BD) is a popular technique to decompose optimization problems (Benders, 1962; Rahmaniani et al., 2017). Contrary to the previous methods, it aims to split an optimization model into two different subsets based on their complicating variables. When the values of those complicating variables are fixed, the resultant problem is less complicated to solve. Thus, by applying BD to our wildfire suppression model, the following is obtained: a master problem, which manages to select what resources can contain the wildfire (complicating variables); and a subproblem, whose goal is to search for feasible solutions according to rest policy while knowing the state of resources in every moment. We have adapted the BD method by adding a set of transformations to the original problem based on the next premise: if it is known when a resource starts to work in the wildfire, we can determine the rest, flight, and work periods. Consequently, to take advantage of this idea, it is necessary to modify the original problem to facilitate separability. This new reformulation will consider how resources must perform rest and travel periods from the start period until the last period (period m). Furthermore, this information is kept in mind to establish when resources must rest, travel, and work while they are in use, that is, until variable eit takes the value 1. Reformulated Problem Transformations can be classified into two groups. First, variables associated with travel times have been divided into three new variables to define traveling periods at the start, rest and end of work (see Table 3). Table 3: Definition of travel variables. Variables Definition trs it ∈ {0,1}It takes the value 1 if resource i∈ I is travelling to go to the wildfire in period t∈ T (travelling associated with start period). trr it ∈ {0,1}It takes the value 1 if resource i∈ I is travelling to perform a rest period or to return to the wildfire in period t∈ T (travelling associated with breaks). tre it ∈ {0,1}It takes the value 1 if resource i∈ I is travelling to leave the wildfire in period t∈ T (travelling associated with end period). Therefore, we know when a resource starts to work in the wildfire and are able to determine the rest, flight, and work periods during all periods. 13
The second transformations are because knowing the start period of a resource can allow one to establish how the resource should act in each period (work, travel, or rest). To accomplish this aim, a new group of variables has been introduced in the model to indicate how a resource must act from its start period until the last one (period ). These variables are denoted by a hat over them: traveling associated with rest periods ( ˆ trr it), resting (ˆrit), ending break (ˆerit) and ending (ˆeit). Auxiliary variables ˆuit,ˆwit and ˆcrit are also defined using the expression of the original model but replacing trit,rit,erit and eit with ˆ trr it,ˆrit,ˆerit and ˆeit, respectively. As a result of applying these transformations, the reformulated problem can be expressed as: minimize (1) subject to (2), (3), (5), (6), (13), (14), (15), (16), (17), (18), (19), ∀i∈ I, t ∈ T , Ai·wit ≤X t0∈T t trs it0(r4) ∀i∈ I, t ∈ T ,X t0∈T t t−T RPi+1 tre it0≥T RPi·eit (r7) ∀i∈ I, t ∈ T ,0≤ˆcrit ≤W Pi(r8) ∀i∈ I, t ∈ T ,ˆrit ≤X t0∈T t+RPi−1 t ˆerit0(r9) ∀i∈ I, t ∈ T :t≥RPi,X t0∈T t t−RPi+1 ˆrit0≥RPi·ˆerit (r10) ∀i∈ I, t ∈ T :t < RPi, CRPi·si1+X t0∈T t ˆrit0≥RPi·ˆerit (r11) ∀i∈ I, t ∈ T ,X t0∈T t+T RPi t−T RPi (ˆrit0+ˆ trr it0)≥X t0∈T t+T RPi t−T RPi ˆrit (r12) ∀i∈ I, t ∈ T ,ˆrit +ˆ trr it ≤ˆuit (r18) Moreover, the following constraints are introduced to establish a relation between the original variables and the new variables: ∀i∈ I, t ∈ T ,uit ≤ˆuit (26) 14
∀i∈ I, t ∈ T ,rit =uit ·ˆrit (27) ∀i∈ I, t ∈ T ,erit =uit ·ˆerit (28) ∀i∈ I, t ∈ T ,wit =uit ·ˆwit ·(1 −max{trs it,tre it})(29) ∀i∈ I, t ∈ T ,trr it =uit ·ˆ trr it (30) ∀i∈ I, t ∈ T ,trit = max{trs it,tre it,trr it}(31) In order to explain the reformulated problem, we illustrate the performed transformations with an example. Example 1. Let us consider a problem instance with a single resource (I= {1}) and 9 periods (T={1, . . . , 9}). Suppose that the resource starts without initial conditions (CW P1=CRP1=CUP1=IT W1=IOW1= 0) and it must travel 1 period from its origin to the wildfire (A1= 1). The resource performance is 1 for all periods (P R1t= 1 for all t∈ T ), and the maximum number of allowed periods without breaks is 4 (W P1= 4). Furthermore, it must perform 1 travel period between rests (T RP1= 1), and 1 rest period (RP1= 1). In this context, suppose that the wildfire has an initial perimeter of 2 km and it grows 0.1 km per period (P ER1= 2 and P ERt= 0.1for all t∈ {2, . . . , 9}). Figure 1 shows the solution (active variables for each period) of the given instance to represent the main idea of the reformulated problem. The figure illustrates how the new variables split the model into two parts to satisfy the rest policy. •The new variables (those over the edges) denoted with the hat represent how the resource must perform the rest periods. Furthermore, it also considers travel due to rest periods. •The original variables (those below the edges) represent how the resource works in the wildfire to contain it, ensuring that the travel and rest periods established by the new variables must be performed. Figure 1 shows that in the first period, the resource starts to be used in the wildfire. At this moment, although variable ˆwindicates that the resource could work, it does not do so because constraints (r4) and (29) force the resource to perform a starting flight before beginning to work. In periods 2, 3, and 7, the 15
0123456789 s=1 ˆw=1 trs=1 s=1 ˆw=1 w=1 ˆw=1 w=1 ˆ trr=1 trr=1 ˆr=1 r=1 ˆ trr=1 trr=1 ˆw=1 w=1 ˆw=1 tre=1 e=1 ˆe=1 ˆ trr=1 Figure 1: Illustration of the relation between original variables and the new ones. resource works because of ˆw= 1, and it does not have to perform starting or ending flights (constraint (29)). In periods 4 and 6, the resource is forced to fly to take a break due to constraint (30). Something similar happens in period 5 with constraints (27) and (28). Finally, the resource contains the wildfire in period 7 with the ending flight occurring in the next period (allowed by the same reasons as the starting flight). It is important to note that the starting period is the same for both cases since the definitions of the auxiliary variables ˆuare defined using s. Table 4 represents the evolution of the wildfire and the performance of the resource over the periods. Note that the resource contains the wildfire in period 7, so from this period, the wildfire perimeter will be 0. In period 8, the resource could keep working since it would satisfy the rest policy ( ˆw= 1), but the wildfire is contained, so the resource must leave it (tre= 1). To simplify the notation, in Table 4, we denote F ireP ertas the perimeter of the wildfire and ResoP ert as the perimeter performed by the resources in each period, i.e., for all t∈ T . F ireP ert:= X t0∈T t P ERt0·yt−1 ResoP ert:= X i∈I,t0∈T t P Rit0·wit0 1 2 3 4 5 6 7 8 9 F ireP ert2.0 2.1 2.2 2.3 2.4 2.5 2.6 0.0 0.0 ResoP ert0.0 1.0 2.0 2.0 2.0 2.0 3.0 3.0 3.0 Table 4: Evolution of the wildfire perimeter and the contention of the resources. As shown, the reformulated problem is more complicated to solve than the original problem since it combines integer variables and nonlinear constraints. However, the purpose is not to solve this problem but to improve its decomposition. The reformulated problem is transformed so that when Benders decomposition is applied, the nonlinear constraints (27)-(30) are linearized by fixing the variables of the original problem. Similar to the notation of a solution in Section 2.2, we define a solution for the original problem as: SOL∗= (s∗ 11, . . . , s∗ nm,tr∗ 11, . . . , tr∗ nm,r∗ 11, . . . , r∗ nm,er∗ 11, . . . , er∗ nm, e∗ 11, . . . , e∗ nm,y∗ 1, . . . , y∗ m,µ∗ 11, . . . , µ∗ gm). 16
In the case of the reformulated problem, a solution is defined similarly by adding values of the new variables at the end of the original solution vector: ˆ SOL∗= (s∗ 11, . . . , s∗ nm,tr∗ 11, . . . , tr∗ nm,r∗ 11, . . . , r∗ nm,er∗ 11, . . . , er∗ nm, e∗ 11, . . . , e∗ nm,y∗ 1, . . . , y∗ m,µ∗ 11, . . . , µ∗ gm,trs∗ 11, . . . , trs∗ nm, tre∗ 11, . . . , tre∗ nm,trr∗ 11, . . . , trr∗ nm, ˆ trr∗ 11, . . . , ˆ trr∗ nm, ˆr∗ 11, . . . , ˆr∗ nm,ˆer∗ 11, . . . , ˆer∗ nm,ˆe∗ 11, . . . , ˆe∗ nm) = (SOL∗, SOL∗ R). To show the equivalence of the original problem and the reformulated problem, the following remark is introduced. Remark 2.1.The auxiliary variables wit are nonnegative for all i∈ I and t∈ T . Proof. This remark can be demonstrated by contradiction. First, from constraint (r18), we know that ˆrit +ˆ trr it ≤ˆuit ⇒uit ·ˆrit +uit ·ˆ trr it ≤uit ·ˆuit ⇒ rit +trr it ≤uit (by definition of wit)⇒rit +trr it ≤wit +rit +trit ⇒ wit ≥trr it −trit. Now, for the sake of contradiction, let us suppose that wit <0. Then, trr it −trit <0⇒trit = 1 and trr it = 0 ⇒ trs it = 1 or tre it = 1 (by equation (29)) ⇒ wit =uit ·ˆwit ·0=0≮0. The following proposition proves the equivalence between both problems, the original problem and the reformulated problem. Proposition 2.2. Let SOL∗be a feasible solution of the original problem; then, there exists an associated feasible solution of the reformulated problem, ˆ SOL∗= (SOL∗, SOL∗ R). Furthermore, let ˆ SOL∗= (SOL∗, SOL∗ R)be a feasible solution of the reformulated problem; then, is a feasible solution of the original problem. Proof. Given a feasible solution of the original problem, SOL∗, it is trivial to prove that it has an associated feasible solution in the reformulated problem by considering each i∈ I,t∈ T ,trs it =tre it =trr it =ˆ trr it =tr∗ it,ˆrit =r∗ it, ˆerit =er∗ it and ˆeit =e∗ it. In addition, we will prove that a feasible solution of the reformulated problem has an associated feasible solution in the original problem. Let us start by proving that constraint (r4) is equivalent to constraint (4). Due to constraint (31), for each i∈ I and for each t∈ T , 17
Ai·wit ≤X t0∈T t trs it0≤X t0∈T t trit0 The equivalence between constraints (r7) and (7) is also trivial. In order to prove the equivalence between constraint (r8) and constraint (8), let us consider the case where IT Wi= 0 and IOWi= 0 (the other case is analogous). Then, for each i∈ I and for each t∈ T , we have crit := X t0∈T t (t+ 1 −t0)·sit0−X t0∈T t (t−t0)·eit0−X t0∈T t rit0−X t0∈T t W Pi·erit0 Now, if period t∗∈ T where eit∗= 1 is considered, for all t≤t∗, crit =X t0∈T t (t+ 1 −t0)·sit0−X t0∈T t rit0−X t0∈T t W Pi·erit0 =X t0∈T t (t+ 1 −t0)·sit0−X t0∈T t uit0 ˆrit0−X t0∈T t W Pi·uit0 ˆerit0 =X t0∈T t (t+ 1 −t0)·sit0−X t0∈T t ˆrit0−X t0∈T t W Pi·ˆerit0 =ˆcrit, where the second equation is because of constraints (27) and (28). The third equality is due to the following: •If ˆrit = 1, by definition of ˆwit and Remark 2.1, it is clear that the variable ˆuit = 1. Otherwise, if ˆrit = 0, the equality is trivial. •If ˆerit = 1, then by constraint (r10) (if t≥RPi) or by constraint (r11) (if t < RPi), it can be deduced that ˆrit = 1. Then, applying the previous item, we know that ˆuit = 1. Otherwise, if ˆerit = 0, the equality is trivial. In the cases where ˆuit = 1, considering the definitions of the auxiliary variables uit and ˆuit, we have that ˆuit =X t0∈T t sit0−X t0∈T t−1 ˆeit0= 1 ⇒X t0∈T t sit0= 1 and therefore2 uit =X t0∈T t sit0−X t0∈T t−1 eit0= 1 −X t0∈T t−1 eit0= 1 proving that if ˆrit = 1 or ˆerit = 1, then uit = 1 for all t≤t∗. 2Note that Pt0∈T t−1eit0=0 for all t≤t∗since eit∗= 1 by assumption. 18
Hence, we have proved that constraint (8) holds for all t≤t∗, 0≤crit =ˆcrit ≤W Pi. Otherwise, if t>t∗, the proof is similar since the expression of crit0takes the same value for periods t0> t∗. This is because uit0= 0 for all t0> t∗, crit =X t0∈T t (t+ 1 −t0)·sit0−X t0∈T t (t−t0)·eit0−X t0∈T t rit0−X t0∈T t W Pi·erit0 =X t0∈T t (t+ 1 −t0)·sit0−(t−t∗)−X t0∈T t uit0·ˆrit0−X t0∈T t W Pi·uit0 ˆerit0 =X t0∈T t (t+ 1 −t0)·sit0−(t−t∗)−X t0∈T t ˆrit0−X t0∈T t W Pi·ˆerit0 =X t0∈T t (t+ 1 −t0)·sit0−(X t0∈T t sit0)·(t−t∗)−X t0∈T t ˆrit0−X t0∈T t W Pi·ˆerit0 =X t0∈T t∗ (t∗+ 1 −t0)·sit0−X t0∈T t∗ ˆrit0−X t0∈T t∗ W Pi·ˆerit0 =ˆcrit∗. In the second equality, we use that eit0= 0 for all t06=t∗and constraints (27) and (28). The third equality can be proven using a procedure similar to that used for the case t≤t∗. The fourth equality is because of the definition of uit and the fact that ˆuit ≥0by constraint (18), which implies that Pt0∈T tsit0= 1. Finally, the fifth equation results from the fact that ˆrit = 0 and ˆerit = 0 for all t>t∗. Hence, we have proven that constraint (8) also holds for all t>t∗, 0≤crit =ˆcrit∗≤W Pi. The proof of equivalence related to constraints (r9)-(r12), is similar to those already proved, but the following considerations are important: 1. The equivalence between constraint (r9) and (9) can be proven by analyzing three different situations: t∈(−∞, t∗−RPi+1],t∈(t∗−RPi+1, t∗] and t∈(t∗,∞). For the proof related to constraints (r10) and (r11) one must distinguish two cases: t∈(−∞, t∗]and t∈(t∗,∞). Finally, for constraint (r12), the proof must be done differentiating between the following cases: t∈(−∞, t∗−T RPi],t∈(t∗−T RPi, t∗+T RPi]and t∈(t∗+T RPit, ∞). 2. Furthermore, for constraints (r9) and (r12) it is necessary to consider the optimality of the solution to demonstrate the cases in which t≤t∗and t≤t∗+T RPi, respectively. Finally, the equivalence between constraints (r18) and (18) is trivial using Remark 2.1: 0≤wit =uit −rit −trit ⇒uit ≥rit +trit. 19
Once the equivalence between the original problem and its reformulation has been demonstrated, we proceed to apply the BD approach to the reformulated problem. Benders Master Problem The master problem seeks the containment of the forest fire without acknowledging the rest periods of the resources: minimize (1) subject to (2), (3), (r4), (5), (6), (r7), (13), (14), (15), (16), (17), (18), (19), (31) ∀(i, t)∈ S∗, m ·sit +X t0∈W− it wit0≤m, (32) ∀(i, t)∈ S∗, m ·sit +X t0∈T R− it trr it0≤m, (33) ∀(i, t)∈ S∗, m ·sit +X t0∈R− it rit0≤m. (34) The variables of the Benders master Problem are the original variables: sit, rit,erit,wit,trs it,trr it and tre it. Moreover, S∗is the set of all the tuples that represent the resources and periods, (i, t)∈ I × T , where resource istarts in period tat some iteration of the algorithm, i.e., S∗:= [ ν∈N S(ν), being S(ν) := {(i, t)∈ I × T :sit(ν)=1}. Using the definition of S∗, cuts (32), (33) and (34) impose where a given (i, t)∈ S∗(iis a resource and tis a starting period) cannot work, rest and travel due to rest periods, respectively. This can be accomplished using the definition of the following sets: W− it :={t0∈ T : ˆwit0(t)=0}, T R− it :={t0∈ T :ˆ trr it0(t)=0}, R− it :={t0∈ T : ˆrit0(t)=0}, where ˆwit0(t),ˆ trr it0(t)and ˆrit0(t)are the values of working, travel due to rest, and rest variables if resource istarts in period t, respectively. 20
It is important to note that the needed information to build cuts (32), (33) and (34) is obtained by solving the Benders subproblem defined below. Benders Subproblem The subproblem seeks to establish a correct policy of breaks considering when resources begin to work. As previously stated, knowing the beginning of the resources allows computing the break periods easily. Thus, the problem is greatly simplified: maximize X i∈I,t∈T ˆwit subject to (r8), (r9), (r10), (r11), (r12), (r18), (26), (27), (28), (29), (30). In the subproblem we employ the previous definition of the variables ˆeit,ˆrit, ˆ trr it and ˆerit, and the expressions defining ˆuit,ˆwit and ˆcrit. Additionally, the variables associated with the master problem (complicating variables) are fixed: sit =s∗ it,rit =r∗ it,erit =er∗ it,wit =w∗ it,trs it =trs∗ it ,trr it =trr∗ it and tre it =tre∗ it . Therefore, as these variables are fixed, constraints (27)-(30) become linear. Remark 2.3.Since the objective function is fully defined in the master problem, in the subproblem, we consider that resources must work all they can. This is a necessary condition for the correct performance of the rest and for wildfire containment. Note that the subproblem only accounts for the feasibility of the solution given by the master problem. Remark 2.4.As the structure of the subproblem shows, it is easy to check that the subproblem can be decomposed trivially by resources. We can do this because the constraints and objective functions do not share information on more than one resource. Thus, instead of solving the Benders subproblem, we can solve for each selected resource i(given from a solution of the Benders master problem). In addition, in order to obtain information on the resolution of these subproblems, we can remove constraints (27)-(30) to obtain information about the configuration of the rest periods. By doing so, the subproblem associated with resource i, which we will name Benders subproblem i, determines the maximum performance of resource iconsidering when it starts and the restrictions on rest period legislation. Benders Algorithm The pseudocode presented in Algorithm 1 defines the steps in the Benders algorithm in the context of our decomposition of the problem. 21
Algorithm 1 Adaptation of the Benders Algorithm 1: ν←0 2: while no new cuts do 3: ν←ν+ 1 4: Solve Benders Master Problem 5: Get S(ν) 6: for (i, t)in S(ν)do 7: if (i, t)not in S∗then 8: Solve Benders Subproblem i 9: Update S∗,W− it ,T R− it and R− it 10: Add the associated cuts given by constraints (32), (33) and (34) It is important to note that three new types of cuts, (32), (33) and (34), are created to improve the convergence of the algorithm. The following remark will discuss the validity of the proposed cuts. Remark 2.5.Constraints (32)-(34) are valid Benders cuts that allow us to give to the Benders master problem information of how a resource must act knowing its start period. For a given i∈ I and ts∈ T where sits= 1, Benders subproblem igives the configuration of work, travel associated with rest, and rest periods of resource i from its start period (ts) until the last period (m). The configuration of work, travel and rest periods guarantees the feasibility of the Benders subproblem i (and consequently, the feasibility of the Benders subproblem). Examining Benders master problem, constraint (32) fixes the working periods to 0according to the solution of Benders subproblem iwith sits= 1. In a similar way, we can impose periods of no travel due to rest and periods of no rest using constraints (33) and (34), respectively. If a solution proposed by the Benders master problem is infeasible for the reformulated problem, then the cuts generated in the next iteration will remove the solution as each Benders subproblem iwill compute a feasible configuration for the resource ito construct the new cuts. Otherwise, if the solution is feasible for the reformulated problem, then it will be feasible for each Benders subproblem iand the associated cuts will admit it. 3. Fixed Activity Problem Inspired by the idea of BD to eliminate the difficulties caused by the computation of the periods in which resources must rest, a new model was developed that guarantees convergence to the global optimum. If the period in which a resource begins to work in the tasks of extinguishing fires and the initial conditions in which this resource starts to work are known, it is simple to compute the work and rest periods while analyzing the pertinent 22
Caunhye, A.M., Nie, X., Pokharel, S., 2012. Optimization models in emergency logistics: A literature review. Socio-Economic Planning Sciences 46, 4–13. Conejo, A., Castillo, E., Mínguez, R., García-Bertrand, R., 2006. Decomposition techniques in mathematical programming: engineering and science applications. Springer Science & Business Media. Dennison, P.E., Brewer, S.C., Arnold, J.D., Moritz, M.A., 2014. Large wildfire trends in the western united states, 1984–2011. Geophysical Research Letters 41, 2928–2933. Dolan, E.D., Moré, J.J., 2002. Benchmarking optimization software with performance profiles. Mathematical Programming 91, 201–213. European Commision, 2019. Forest fires in Europe, Middle East and North Africa 2018. Technical Report. JRC Technical Reports. Fisher, M.L., 2004. The lagrangian relaxation method for solving integer programming problems. Management science 50, 1861–1871. Gleixner, A., Bastubbe, M., Eifler, L., Gally, T., Gamrath, G., Gottwald, R.L., Hendel, G., Hojny, C., Koch, T., Lübbecke, M.E., Maher, S.J., Miltenberger, M., Müller, B., Pfetsch, M.E., Puchert, C., Rehfeldt, D., Schlösser, F., Schubert, C., Serrano, F., Shinano, Y., Viernickel, J.M., Walter, M., Wegscheider, F., Witt, J.T., Witzig, J., 2018. The SCIP Optimization Suite 6.0. Technical Report. Optimization Online. URL: http://www.optimization-online. org/DB_HTML/2018/07/6692.html. Gurobi Optimization, L., 2020. Gurobi optimizer reference manual. URL: http: //www.gurobi.com. Hamdi, A., Mishra, S.K., 2011. Decomposition methods based on augmented lagrangians: a survey, in: Topics in nonconvex optimization. Springer, pp. 175–203. Lee, W., 2006. A stochastic mixed integer programming approach to wildfire management systems. Texas A&M University. Martell, D.L., 2007. Forest fire management, in: Handbook of operations research in natural resources. Springer, pp. 489–509. Martins, I., Alvelos, F., Constantino, M., 2012. A branch-and-price approach for harvest scheduling subject to maximum area restrictions. Computational Optimization and Applications 51, 363–385. Mercier, A., Cordeau, J.F., Soumis, F., 2005. A computational study of benders decomposition for the integrated aircraft routing and crew scheduling problem. Computers & Operations Research 32, 1451–1476. Miller, C., Ager, A.A., 2013. A review of recent advances in risk analysis for wildfire management. International Journal of Wildland Fire 22, 1–14. 29
Minas, J.P., Hearne, J.W., Handmer, J.W., 2012. A review of operations research methods applicable to wildfire management. International Journal of Wildland Fire 21, 189–196. Nishi, T., Ando, M., Konishi, M., 2005. Distributed route planning for multiple mobile robots using an augmented lagrangian decomposition and coordination technique. IEEE Transactions on Robotics 21, 1191–1200. Papadakos, N., 2009. Integrated airline scheduling. Computers & Operations Research 36, 176–195. Rahmaniani, R., Crainic, T.G., Gendreau, M., Rei, W., 2017. The benders decomposition algorithm: A literature review. European Journal of Operational Research 259, 801–817. doi:10.1016/j.ejor.2016.12.005. Rios, J., Ross, K., 2010. Massively parallel dantzig-wolfe decomposition applied to traffic flow scheduling. Journal of Aerospace Computing, Information, and Communication 7, 32–45. Rodríguez-Veiga, J., Ginzo-Villamayor, M.J., Casas-Méndez, B., 2018a. An integer linear programming model to select and temporally allocate resources for fighting forest fires. Forests 9, 583. Rodríguez-Veiga, J., Gómez-Costa, I., Ginzo-Villamayor, M.J., Casas-Méndez, B., Sáiz-Díaz, J.L., 2018b. Assignment problems in wildfire suppression: Models for optimization of aerial resource logistics. Forest Science 64, 504–514. Romanski, J., Hentenryck, P.V., 2016. Benders decomposition for large-scale prescriptive evacuations. Van Wassenhove, L.N., Pedraza Martinez, A.J., 2012. Using or to adapt supply chain management best practices to humanitarian logistics. International Transactions in Operational Research 19, 307–322. Vanderbeck, F., Savelsbergh, M.W., 2006. A generic view of dantzig–wolfe decomposition in mixed integer programming. Operations Research Letters 34, 296–306. Vanderbeck, F., Wolsey, L.A., 1996. An exact algorithm for IP column generation. Operations Research Letters 19, 151–159. Yang, Z., Guo, L., Yang, Z., 2019. Emergency logistics for wildfire suppression based on forecasted disaster evolution. Annals of Operations Research 283, 917–937. Zhou, S., Erdogan, A., 2019. A spatial optimization model for resource allocation for wildfire suppression and resident evacuation. Computers & Industrial Engineering 138, 106101. Úbeda, X., Sarricolea, P., 2016. Wildfires in chile: A review. Global and Planetary Change 146, 152–161. 30
A. Dimensions of the instances Brigades Aircraft Machines Periods Variables Constraints 2 2 2 10 340 521 2 2 2 15 510 766 2 2 4 10 440 671 2 2 4 15 660 986 2 4 2 10 440 671 2 4 2 15 660 986 2 4 4 10 540 821 2 4 4 15 810 1206 4 2 2 10 440 671 4 2 2 15 660 986 4 2 4 10 540 821 4 2 4 15 810 1206 4 4 2 10 540 821 4 4 2 15 810 1206 4 4 4 10 640 971 4 4 4 15 960 1426 Table A.5: Dimensions of the problems. 31
Brigades Aircraft Machines Periods Variables Constraints 5 5 5 10 790 1196 5 5 5 20 1580 2316 5 5 5 30 2370 3436 5 5 5 40 3160 4556 5 5 5 50 3950 5676 5 5 5 60 4740 6796 5 5 20 10 1540 2321 5 5 20 20 3080 4491 5 5 20 30 4620 6661 5 5 20 40 6160 8831 5 5 20 50 7700 11001 5 5 20 60 9240 13171 5 20 5 10 1540 2321 5 20 5 20 3080 4491 5 20 5 30 4620 6661 5 20 5 40 6160 8831 5 20 5 50 7700 11001 5 20 5 60 9240 13171 5 20 20 10 2290 3446 5 20 20 20 4580 6666 5 20 20 30 6870 9886 5 20 20 40 9160 13106 5 20 20 50 11450 16326 5 20 20 60 13740 19546 20 5 5 10 1540 2321 20 5 5 20 3080 4491 20 5 5 30 4620 6661 20 5 5 40 6160 8831 20 5 5 50 7700 11001 20 5 5 60 9240 13171 20 5 20 10 2290 3446 20 5 20 20 4580 6666 20 5 20 30 6870 9886 20 5 20 40 9160 13106 20 5 20 50 11450 16326 20 5 20 60 13740 19546 20 20 5 10 2290 3446 20 20 5 20 4580 6666 20 20 5 30 6870 9886 20 20 5 40 9160 13106 20 20 5 50 11450 16326 20 20 5 60 13740 19546 20 20 20 10 3040 4571 20 20 20 20 6080 8841 20 20 20 30 9120 13111 20 20 20 40 12160 17381 20 20 20 50 15200 21651 20 20 20 60 18240 25921 Table A.6: Dimensions of the problems in large cases. 32
B. Auxiliary Figures 20-20-20-10 20-20-20-20 20-20-20-30 20-20-20-40 20-20-20-50 20-20-20-60 20-20-5-10 20-20-5-20 20-20-5-30 20-20-5-40 20-20-5-50 20-20-5-60 20-5-20-10 20-5-20-20 20-5-20-30 20-5-20-40 20-5-20-50 20-5-20-60 20-5-5-10 20-5-5-20 20-5-5-30 20-5-5-40 20-5-5-50 20-5-5-60 5-20-20-10 5-20-20-20 5-20-20-30 5-20-20-40 5-20-20-50 5-20-20-60 5-20-5-10 5-20-5-20 5-20-5-30 5-20-5-40 5-20-5-50 5-20-5-60 5-5-20-10 5-5-20-20 5-5-20-30 5-5-20-40 5-5-20-50 5-5-20-60 5-5-5-10 5-5-5-20 5-5-5-30 5-5-5-40 5-5-5-50 5-5-5-60 0 100 200 300 400 500 600 FA OR Elapsed Time Figure B.6: Boxplot of the computational time of the algorithm only considering the solved instances. X axis represent the size of the instances, representing the numbers separated by hyphens the number of brigades, aircraft, machines and periods, respectively. 33
CAPÍTULO 5. TRABAJOS PUBLICADOS Y SOMETIDOS A REVISIÓN 5.4. Assignment problems in wildfire suppression: models for optimization of aerial resource logistics Referencia del artículo J. Rodríguez-Veiga, I. Gómez-Costa, M. J. Ginzo-Villamayor y col., «Assignment problems in wildfire suppression: Models for optimization of aerial resource logistics,» Forest Science, vol. 64, n.o5, págs. 504-514, 2018 Filiación autores: Jorge Rodríguez-Veiga1, Iván Gómez-Costa2, María José Ginzo-Villamayor3, Balbina Casas-Méndez34, José Luís Sáiz-Díaz5. Contribución: Modelado, implementación, simulaciones y redacción del artículo. Factor de impacto: 1.058. Categoría: Forestry-scie. Posición relativa: Nº44 de un total de 67 (Q3). Citas scopus: 3. Citas Google: 8. ISSN: 1938-3738. Enlace: https://academic.oup.com/forestscience/ 1Departamento de Estadística, Análisis Matemático y Optimización, Universidad de Santiago de Compostela, Santiago de Compostela, España. 2INDRA, Technology and Consulting, A Coruña, España. 3Grupo de Investigación Modestya, Departamento de Estadística, Análisis Matemático y Optimización, Universidad de Santiago de Compostela, Santiago de Compostela, España. 4Facultad de Matemáticas, Campus Vida s/n 15782 Santiago de Compostela, España. 5INAER, Babcock Internacional España, Alicante, España. 87
Jorge Rodríguez Veiga 5.5. Wildfire resources management: a decision support tool created with R to solve optimization models in logistics for fighting forest fires Referencia del capítulo J. Rodríguez Veiga, M. J. Ginzo Villamayor y B. V. Casas Méndez, «Wildfire resources management: a decision support tool created with R to solve optimization models in logistics for fighting forest fires,» en Progress in Industrial Mathematics: Success Stories, M. Cruz, C. Parés y P. Quintela, eds., vol. 5, Springer International Publishing, 2021, págs. 227-246 Filiación autores: Jorge Rodríguez-Veiga12345, María José Ginzo-Villamayor12345, Balbina CasasMéndez12345. Contribución: Modelado, implementación, simulaciones y revisión de la redacción del capítulo. Citas scopus: 0. Citas Google: 0. ISSN: 2662-7183. Enlace: https://www.springerprofessional.de 1Departamento de Estadística, Análisis Matemático y Optimización 2Instituto de Matemáticas (IMAT), España. 3Instituto Tecnológico de Matemática Industrial (ITMATI), España. 4Grupo de Modelos de Optimización, Decisión, ESTadística Y Aplicaciones (MODESTYA), España. 5Universidad de Santiago de Compostela, España. 88
Capítulo 6 Conclusiones y trabajo futuro El trabajo realizado en la tesis contempla el cálculo de soluciones en problemas de investigación operativa. En concreto se estudia el cálculo de soluciones en dos ámbitos, la teoría de juegos, en juegos de mayoría ponderada y mayoría ponderada múltiple con configuración de coaliciones, y la optimización, en la gestión de recursos en la contención de incendios forestales. En la Sección 6.1 se describen las conclusiones y trabajo futuro respecto al cálculo de soluciones en juegos de mayoría ponderada múltiple. En la Sección 6.2 se describen las conclusiones y trabajo futuro respecto al cálculo de soluciones para una gestión eficiente de los recursos que participan en la contención de un incendio forestal. En ambas secciones se describen las conclusiones del trabajo realizado y las posibles líneas de trabajos para continuar las investigaciones realizadas. Cabe resaltar que todos los trabajos realizados fueron acompañados de la programación de los modelos matemáticos o algoritmos allí descritos. Además, se pretendió profundizar en los detalles matemáticos aplicando y demostrando resultados teóricos sobre los problemas planteados. En el caso de los índices de poder, mediante el empleo de las funciones generatrices, y en el caso de la gestión de recursos en la contención de incendios forestales, mediante la aplicación de técnicas de descomposición y reformulaciones del problema para mejorar la eficiencia en la resolución del problema SARC. 89
Jorge Rodríguez Veiga Por último, señalar que en función de las necesidades que fueron surgiendo a lo largo de la tesis, se optó por el empleo de diferentes lenguajes de modelado algebraico y de programación. 6.1. Índices de poder Dado que Albizuri y Aurrekoetxea [3] y Albizuri, Aurrekoetxea y Zarzuelo [2] introdujeron el índice de Banzhaf-Coleman generalizado y el índice con configuración, respectivamente, para juegos simples con configuraciones de coalición, como complemento a la definición y el estudio de sus propiedades realizados en estos documentos, una pregunta abierta era el cálculo eficiente de estos índices mediante las llamadas funciones generatrices, en línea con estudios similares relacionados. En la Sección 5.1 se presenta el cálculo del índice con configuración y Banzhaf-Coleman generalizado mediante el empleo de funciones generatrices. Se presenta el cálculo de los índices demostrando matemáticamente su idoneidad, completando el artículo con la implementación del algoritmo utilizando una herramienta de software libre como R[65] y la presentación de ejemplos de la vida real que muestra el alcance del modelo considerado y los algoritmos introducidos. Por último se realiza una extensión de los algoritmos a la clase más amplia de juegos de mayoría ponderada múltiple. Creemos que podría ser de interés extender estas técnicas a clases más amplias de juegos como los juegos simples monótonos, comenzando por los que no contienen restricciones en la comunicación. Para esta tarea, puede ser de ayuda el trabajo de Freixas y Puente [37]. En este trabajo, se establece explícitamente lo siguiente: “This question is solved in this paper, by means of a constructive proof, which provides a way to represent any complete simple game with minimum as the intersection of a number mof weighted majority games, where mcoincides with its dimension”. Sin embargo, para los juegos de utilidad general transferibles, existe otra técnica de cálculo habitual que podría explorarse en el contexto actual y se basa en las llamadas extensiones multilineales (ver, por ejemplo, los artículos de Owen [62] y [61]). Por otra parte, también cabe señalar que el modelo de juegos cooperativos con configuraciones 90
CAPÍTULO 6. CONCLUSIONES Y TRABAJO FUTURO de coalición estudiado en este trabajo, está siendo de gran interés tal y como recogen los trabajos de Albizuri y Vidal-Puga [4] o Andjiga y Courtin [9], entre otros. Por último, también cabría considerar nuevos tipos de restricciones en la comunicación de los jugadores. Por ejemplo, la estructura de configuración de coaliciones, podría verse enriquecida con la inclusión de jugadores incompatibles [57]. Por todo ello, se cree interesante explorar juegos con nuevas restricciones en la comunicación, empleando para el cálculo de sus índices los algoritmos introducidos en esta tesis u otras herramientas de análisis combinatorio como las propuestas en el trabajo de Neto [56]. También cabe destacar, como trabajo futuro, la necesidad de realizar un estudio de la complejidad de los algoritmos definidos en la Sección 5.1. Es importante saber la complejidad del algoritmo para poder realizar una comparativa entre distintos métodos de cálculo de soluciones. Además, sería de gran utilidad un análisis de los tiempos de computación entre distintos algoritmos de cálculo de soluciones en el contexto establecido. En muchas ocasiones, debido al gran número de jugadores involucrados en procesos de votación reales, estos cálculos se hacen de forma muy ineficiente. Es por ello que realizar un estudio de la eficiencia de los métodos podría aportar gran valor. 6.2. Gestión de recursos en la contención de incendios forestales En las Secciones 5.2 y 5.4 se describen tres modelos de investigación operativa para la gestión de recursos en incendios forestales. Los modelos propuestos surgen de la colaboración con la empresa Babcock España. El modelo propuesto para la resolución del problema SARC extiende el propuesto por Donovan y Rideout [29]. Además, pese a que los problemas AARV yAAPR surgen de forma habitual en la gestión de incendios forestales, no se encontró literatura donde se les diese solución. De los problemas anteriormente mencionados, en esta tesis se establece la relación que existe entre los mismos. Así, se genera conocimiento susceptible de ser transferido a las empresas encargadas de la extinción de 91
Jorge Rodríguez Veiga [29] G. H. Donovan y D. B. Rideout, «An integer programming model to optimize resource allocation for wildfire containment,» Forest Science, vol. 49, n.o2, págs. 331-335, 2003. [30] I. Dragan, T. Driessen e Y. Funaki, «Collinearity between the Shapley value and the egalitarian division rules for cooperative games,» Operations-Research-Spektrum, vol. 18, n.o2, págs. 97-105, 1996. [31] P. Dubey, «On the uniqueness of the Shapley value,» International Journal of Game Theory, vol. 4, n.o3, págs. 131-139, 1975. [32] P. Dubey y L. S. Shapley, «Mathematical properties of the Banzhaf power index,» Mathematics of Operations Research, vol. 4, n.o2, págs. 99-131, 1979. [33] V. Feltkamp, «Alternative axiomatic characterizations of the Shapley and Banzhaf values,» International Journal of Game Theory, vol. 24, n.o2, págs. 179-186, 1995. [34] M. A. Finney, FARSITE: Fire Area Simulator - Model Development and Evaluation. United States Department of Agriculture, Forest Service, Rocky Mountain Research Station, 1998, vol. 3. [35] C. A. Floudas, Nonlinear and mixed-integer optimization: fundamentals and applications. Oxford University Press, 1995. [36] R. Fourer, D. M. Gay y B. W. Kernighan, «A modeling language for mathematical programming,» Management Science, vol. 36, n.o5, págs. 519-554, 1990. [37] J. Freixas y M. A. Puente, «Dimension of complete simple games with minimum,» European Journal of Operational Research, vol. 188, n.o2, págs. 555-568, 2008. [38] A. M. Geoffrion, «Generalized Benders decomposition,» Journal of Optimization Theory and Applications, vol. 10, n.o4, págs. 237-260, 1972. 98
CAPÍTULO 6. CONCLUSIONES Y TRABAJO FUTURO [39] D. B. Gillies, «Solutions to general non-zero-sum games,» en Contributions to the Theory of Games, A. W. Tucker y R. D. Luce, eds., vol. 4, Princeton University Press, 1959, págs. 47-85. [40] Gobierno de España. Ministerio de Agricultura, Pesca y Alimentación, Los Incendios Forestales en España: 1 enero – 31 diciembre 2018 Avance Informativo,https : / / www . mapa . gob.es/es/desarrollo-rural/estadisticas/iiff_2018_ tcm30-507741.pdf, Accessed:14 de junio de 2021, 2019. [41] J. K. Gorte y R. W. Gorte, «Application of economic techniques to fire management-a status review and evaluation,» Gen. Tech. Rep. INT-GTR-53. Ogden, UT: US Department of Agriculture, Forest Service, Intermountain Research Station. 26 p., vol. 53, 1979. [42] L. Gurobi Optimization, Gurobi Optimizer Reference Manual, Available at http://www.gurobi.com, 2020. [43] A. Hamdi y S. K. Mishra, «Decomposition methods based on augmented Lagrangians: a survey,» en Topics in nonconvex optimization, Springer, 2011, págs. 175-203. [44] R. Headley, «Fire Suppression, District 5,» USDA-Forest Service. 57 p., págs. 1-57, 1916. [45] L. S. Lasdon, Optimization theory for large systems. Courier Corporation, 2002. [46] E. L. Lawler y D. E. Wood, «Branch-and-bound methods: A survey,» Operations Research, vol. 14, n.o4, págs. 699-719, 1966. [47] C. Lemaréchal, «Lagrangian relaxation,» en Computational combinatorial optimization, Springer, 2001, págs. 112-156. [48] W. F. Lucas, «Measuring power in weighted voting systems,» en Political and related models, Springer, 1983, págs. 183-238. [49] D. G. Luenberger, Y. Ye y col., Linear and nonlinear programming. Springer, 1984, vol. 2. 99
Jorge Rodríguez Veiga [50] T. L. Magnanti, P. Mireault y R. T. Wong, «Tailoring Benders decomposition for uncapacitated network design,» en Netflow at Pisa, G. Gallo y C. Sandi, eds., Springer, 1986, págs. 112-154. [51] T. L. Magnanti y R. T. Wong, «Accelerating Benders decomposition: Algorithmic enhancement and model selection criteria,» Operations Research, vol. 29, n.o3, págs. 464-484, 1981. [52] D. L. Martell, «A review of recent forest and wildland fire management decision support systems research,» Current Forestry Reports, vol. 1, n.o2, págs. 128-137, 2015. [53] T. Matsui e Y. Matsui, «A survey of algorithms for calculating power indices of weighted majority games,» Journal of the Operations Research Society of Japan, vol. 43, n.o1, págs. 71-86, 2000. [54] J. F. McCloskey, «OR Forum–the beginnings of operations research: 1934–1941,» Operations Research, vol. 35, n.o1, págs. 143-152, 1987. [55] C. Miller y A. A. Ager, «A review of recent advances in risk analysis for wildfire management,» International Journal of Wildland Fire, vol. 22, n.o1, págs. 1-14, 2013. [56] A. F. Neto, «Generating functions of weighted voting games, MacMahon’s partition analysis, and Clifford algebras,» Mathematics of Operations Research, vol. 44, n.o1, págs. 74-101, 2019. [57] A. F. Neto y C. R. Fonseca, «An approach via generating functions to compute power indices of multiple weighted voting games with incompatible players,» Annals of Operations Research, vol. 279, n.o1, págs. 221-249, 2019. [58] J. v. Neumann, «Zur theorie der gesellschaftsspiele,» Mathematische Annalen, vol. 100, n.o1, págs. 295-320, 1928. 100
CAPÍTULO 6. CONCLUSIONES Y TRABAJO FUTURO [59] W. Ongsakul y N. Petcharaks, «Unit commitment by enhanced adaptive Lagrangian relaxation,» IEEE Transactions on Power Systems, vol. 19, n.o1, págs. 620-628, 2004. [60] G. Owen, «Modification of the Banzhaf-Coleman index for games with a priori unions,» en Power, voting, and voting power, M. J. Holler, ed., Springer, 1981, págs. 232-238. [61] G. Owen, «Multilinear extensions and the Banzhaf value,» Naval Research Logistics Quarterly, vol. 22, n.o4, págs. 741-750, 1975. [62] G. Owen, «Multilinear extensions of games,» Management Science, vol. 18, n.o5-part-2, págs. 64-79, 1972. [63] G. Owen, «Values of games with a priori unions,» en Mathematical Economics and Game Theory, O. M. R. Henn, ed., Springer, 1977, págs. 76-88. [64] Python Software Foundation, Python Language Reference, Available at http://www.python.org, 2020. [65] R Core Team, R: A Language and Environment for Statistical Computing, Available at https://www.R - project.org/, 2020. [66] J. Rajgopal, «Principles and applications of operations research,» en Maynard’s Industrial Engineering Handbook, K. B. Zandin, ed., McGraw-Hill Education, 2004, págs. 11-27. [67] K. Ríbnikov y K. Medkov, Análisis combinatorio: problemas y ejercicios. Mir, 1988. [68] J. Rodríguez Veiga, M. J. Ginzo Villamayor y B. V. Casas Méndez, «Wildfire resources management: a decision support tool created with R to solve optimization models in logistics for fighting forest fires,» en Progress in Industrial Mathematics: Success Stories, M. Cruz, C. Parés y P. Quintela, eds., vol. 5, Springer International Publishing, 2021, págs. 227-246. [69] J. Rodríguez-Veiga, ROMO, Available at https://github. com/jorgerodriguezveiga/romo, 2020. 101
Jorge Rodríguez Veiga [70] J. Rodríguez-Veiga, M. J. Ginzo-Villamayor y B. CasasMéndez, «An integer linear programming model to select and temporally allocate resources for fighting forest fires,» Forests, vol. 9, n.o10, 583, págs. 1-18, 2018. [71] J. Rodríguez-Veiga, I. Gómez-Costa, M. J. Ginzo-Villamayor, B. Casas-Méndez y J. L. Sáiz-Díaz, «Assignment problems in wildfire suppression: Models for optimization of aerial resource logistics,» Forest Science, vol. 64, n.o5, págs. 504-514, 2018. [72] J. Rodríguez-Veiga, G. I. Novoa-Flores y B. Casas-Méndez, «Implementing generating functions to obtain power indices with coalition configuration,» Discrete Applied Mathematics, vol. 214, págs. 1-15, 2016. [73] J. Rodríguez-Veiga, D. Rodríguez-Penas, Á. M. GonzálezRueda y M. J. Ginzo-Villamayor, «Application of decomposition techniques in a wildfire suppresion optimization model,» inf. téc., 2021. [74] G. K. Saharidis, M. Minoux y M. G. Ierapetritou, «Accelerating Benders method using covering cut bundle generation,» International Transactions in Operational Research, vol. 17, n.o2, págs. 221-237, 2010. [75] S. Sethi y G. Sorger, «A theory of rolling horizon decision making,» Annals of Operations Research, vol. 29, n.o1, págs. 387-415, 1991. [76] L. S. Shapley, «A value for n-person games,» en Contributions to the Theory of Games, A. W. Tucker y H. Kuhn, eds., vol. 2, Princeton University Press, 1953, págs. 307-317. [77] L. S. Shapley y M. Shubik, «A method for evaluating the distribution of power in a committee system,» American Political Science Review, vol. 48, n.o3, págs. 787-792, 1954. 102
CAPÍTULO 6. CONCLUSIONES Y TRABAJO FUTURO [78] C. Sievert, C. Parmer, T. Hocking, S. Chamberlain, K. Ram, M. Corvellec y P. Despouy, plotly: Create Interactive Web Graphics via ’plotly.js’, R package version 4.7.1, 2017. dirección: https://CRAN.R-project.org/package=plotly. [79] Spanish Ministry of Development, Operational Circular 16B,http : / / www . aecaweb . com / informes / documentos / INFORMES_Y_ESTUDIOS/circular_operativa_16_b.doc. Accessed: 14 de junio de 2021, 1995. [80] W. N. Sparhawk, The use of liability ratings in planning forest fire protection. National Emergency Training Center, United States of America, 1925. [81] F. Vanderbeck y L. A. Wolsey, «An exact algorithm for IP column generation,» Operations Research Letters, vol. 19, n.o4, págs. 151-159, 1996. [82] F. Vanderbeck, «On Dantzig-Wolfe decomposition in integer programming and ways to perform branching in a branchand-price algorithm,» Operations Research, vol. 48, n.o1, págs. 111-128, 2000. [83] J. Von Neumann y O. Morgenstern, Theory of games and economic behavior. Princeton University Press, 1944. [84] E. Winter, «The consistency and potential for values of games with coalition structure,» Games and Economic Behavior, vol. 4, n.o1, págs. 132-144, 1992. [85] X. Zhao, P. B. Luh y J. Wang, «Surrogate gradient algorithm for Lagrangian relaxation,» Journal of Optimization Theory and Applications, vol. 100, n.o3, págs. 699-712, 1999. 103
Apéndice A Técnicas de descomposición En el presente capítulo se describen las ideas principales de las técnicas de descomposición empleadas en el trabajo [73], presentado en la Sección 5.3. Las técnicas de descomposición parten de dividir el problema original en problemas de menor dimensión, con el fin de que la resolución del problema de partida sea más sencilla. Se puede decir que este tipo de técnicas siguen la filosofía de “divide y vencerás” y para realizar las descomposiciones, se apoyan en una estructura concreta de los problemas. Los problemas con los que se va a trabajar se caracterizan por tener dos estructuras bien diferenciadas. La primera posee un subconjunto de restricciones, llamadas restricciones complicantes, que en caso de ser omitidas permiten su descomposición en problemas de menor complejidad. La segunda se caracteriza por tener un subconjunto de variables, llamadas variables complicantes, que al ser fijadas a un valor facilitan la resolución del problema resultante. Por ejemplo, son complicantes las variables enteras en un problema MILP, variables que provocan no linealidades en restricciones y que en caso de ser fijadas a un valor proporcionan un problema lineal, o variables que permiten identificar algún tipo de estructura en el problema que permita separarlo en subproblemas más sencillos. Para la explicación de los distintos algoritmos de descomposición será 105
Jorge Rodríguez Veiga necesario el empleo de algunos resultados clásicos de la programación lineal. Teorema A.1 (Carathéodory).Sea X={x∈Rn:A·x≤b, x≥0}un conjunto (poliedro convexo) no vacío, con A∈Rm×n,b∈Rm. Entonces el conjunto de sus puntos extremos es no vacío y finito, {x1, . . . , xp}= {xi∈Rn:i∈A}. Además, su conjunto de direcciones extremas es vacío si y solo si Xes acotado. Si Xes no acotado, entonces el conjunto de direcciones extremas, {d1, . . . , dl}={dj∈Rn:j∈D}, es no vacío y finito. Asimismo, ¯x∈Xsi y solo si puede ser representado como combinación convexa de sus puntos extremos más una combinación lineal no negativa de sus direcciones extremas, ¯x=X i∈A λi·xi+X j∈D µj·dj X i∈A λi= 1 λi≥0,∀i∈A µj≥0,∀j∈D. Teorema A.2 (Fundamental de la Programación Lineal).Si un problema de programación lineal tiene una solución óptima acotada, entonces existe un punto extremo de su región factible que es óptimo. Teorema A.3 (Dualidad Débil).Dado un problema de programación lineal (problema primal) y su problema dual, y un par de soluciones factibles, una de cada problema, el valor del objetivo del problema de minimización es mayor o igual que el valor del objetivo del problema de maximización, cuando se evalúan en las soluciones consideradas. Teorema A.4 (Dualidad Fuerte).Dado un problema de programación lineal (problema primal) y su problema dual, si uno tiene solución óptima finita, entonces, el otro también la tiene, y los valores óptimos de sus funciones objetivo coinciden. Teniendo en cuenta el Teorema de Dualidad Débil A.3, se obtiene el siguiente corolario. Corolario A.5. Dado un problema de programación lineal (problema primal) y su problema dual, si la función objetivo de uno de ellos es no acotada, entonces el otro no tiene soluciones factibles. 106
APÉNDICE A. TÉCNICAS DE DESCOMPOSICIÓN Los siguientes apartados de este apéndice describen tres técnicas de descomposición empleadas para la resolución de problemas con restricciones complicantes o variables complicantes. En las secciones A.1 y A.2 se describen la técnica de descomposición lagrangiana y la técnica de descomposición de Dantzig-Wolfe, respectivamente. Ambas se emplean sobre problemas con restricciones complicantes, permitiendo manejarlas para obtener una descomposición que simplifique la resolución del problema original. En la Sección A.3 se describe la técnica de descomposición de Benders, empleada para facilitar la resolución de problemas con una estructura con variables complicantes. En estas secciones se describirán, sin profundizar, los algoritmos que permiten aplicar cada una de las técnicas de descomposición. En caso de querer ahondar en los fundamentos teóricos, se facilitarán referencias útiles para cada una de las técnicas. A.1. Descomposición lagrangiana Consideremos el siguiente problema de programación lineal con restricciones complicantes. m´ın x1∈Rn1,x2∈Rn2c| 1·x1+c| 2·x2(A.1) s.a: A1·x1≤a1(A.2) A2·x2≤a2(A.3) B1·x1+B2·x2≤b(A.4) donde para cada subconjunto de variables x1yx2,c1∈Rn1yc2∈Rn2 son sus vectores de coeficientes en la función objetivo y A1∈Rm1×n1y A2∈Rm2×n2sus matrices de coeficientes en las restricciones en las que únicamente están presentes variables de uno de los dos grupos. Además, se representa por a1∈Rm1ya2∈Rm2a los vectores de términos independientes de las restricciones (A.2) y (A.3), respectivamente. Las restricciones (A.4) forman el subconjunto de restricciones complicantes, siendo B1∈RmB×n1yB2∈RmB×n2sus matrices de coeficientes para las variables x1yx2, respectivamente, y b∈RmBsu vector de términos independientes. 107
Jorge Rodríguez Veiga situaciones: 1. No tiene solución factible, entonces el problema original tampoco. 2. Tiene óptimo finito, entonces estamos ante un punto extremo de la región factible X. 3. No tiene óptimo finito, entonces obtendremos una dirección extrema djde la región factible Xresolviendo el siguiente subproblema cónico, m´ın d∈Rnc|·d−ω|·E·d s.a: B·d= 0 ~ 1|·d≤1 d≥0, (A.11) donde ~ 1 = (1,...,1) ∈Rn. Para establecer el criterio de parada del algoritmo, cabe destacar que el objetivo en una solución del problema maestro (A.9) es una cota superior del óptimo del problema original (A.7). Para determinar una cota inferior, dada una solución factible del problema original, ¯x, también factible en el subproblema (A.10), se verifica la siguiente desigualdad: c|·¯x−ω|·E¯x≥ZS, lo que implica que, c|·¯x≥ZS+ω|·E¯x=ZS+ω|·e, donde la igualdad anterior se tiene por ser ¯xuna solución factible del problema original y por tanto verificar que E·¯x=e. Al ser c|·¯x≥ZS+ω|·e para cualquier solución factible del problema original, ZS+ω|·eserá una cota inferior del objetivo óptimo. Una vez presentados el problema maestro y el subproblema así como las relaciones entre ellos, pasamos a describir el algoritmo de descomposición de Dantzig-Wolfe. El algoritmo comienza mediante la inicialización del conjunto de puntos y direcciones extremas así como de las cotas inferior y superior para el valor del objetivo. A continuación se resuelve 114
APÉNDICE A. TÉCNICAS DE DESCOMPOSICIÓN el problema maestro (A.9). Se obtienen su solución primal (λ, µ)y su solución dual (ω, α). Se actualiza la cota superior con el valor de la función objetivo. Posteriormente, se considera el subproblema (A.10) y se obtiene un punto extremo (xp) o una dirección extrema (dl) de X. En el primer caso se actualiza la cota inferior si ZS+w|·ees mayor que la cota inferior actual. Si las cotas inferior y superior coinciden el algoritmo finaliza. En otro caso, se incorpora al problema maestro el punto o dirección extrema obtenida y se repiten nuevamente los pasos descritos. El Algoritmo 2 muestra el pseudocódigo de la descomposición de DantzigWolfe. Una explicación más detallada del algoritmo se puede ver en [13] y [49], en el manual orientado a la aplicación [17], y en el libro [27], entre otros. Es importante destacar que el algoritmo requiere partir de un subconjunto de puntos extremos y direcciones extremas (este último podría ser vacío). Para inicializar los subconjuntos existen diversos métodos. Una de ellos consiste en realizar perturbaciones de la función objetivo de los subproblemas (A.10) buscando obtener distintos puntos extremos de su región factible. Otra opción es la resolución del algoritmo de DantzigWolfe pero sobre una reformulación del problema maestro (A.9) y del subproblema (A.10). El problema maestro se reformula incorporando variables auxiliares en las restricciones y omitiendo la función objetivo original y minimizando en esta una penalización por el empleo de las variables auxiliares. El subproblema se reformula omitiendo en este el término c|·x. La descripción detallada de este algoritmo se puede ver en [13]. Además el algoritmo se puede extender para descomponer el subproblema (A.10) en varios subproblemas, siempre y cuando la estructura lo permita. El algoritmo se basa en el procedimiento expuesto, pero en este caso habrá que aplicar el Teorema de Carathéodory A.1 sobre la región factible de cada subproblema. Detalles de la extensión del algoritmo se pueden ver en [13]. Por último, se menciona que, para problemas MILP, la técnica se conoce como branch-and-price [12], [82], y aúna las metodologías de Dantzig-Wolfe y branch-and-bound [46]. 115
Jorge Rodríguez Veiga Algoritmo 2 Descomposición de Dantzig-Wolfe Input: Inicializar LB ← −∞,UB ← ∞ y los subconjuntos de puntos extremos y direcciones extremas ARyDR, respectivamente. 1: while LB < UB do 2: Resolver el problema maestro (A.9). 3: Obtener el valor de su función objetivo, Z, las variables primales (λ, µ)y las variables duales (ω, α). 4: Actualizar UB ←Z. 5: Resolver el subproblema (A.10). 6: Obtener el valor de su función objetivo ZS. 7: if ZS>−∞ then 8: Obtener la variable primal x. 9: Añadir el punto extremo xp←x(con p=|AR|+ 1) al conjunto de puntos extremos AR←AR∪ {p}. 10: Actualizar LB ←m´ax{LB, ZS+ω|·e}. 11: else 12: Resolver el problema (A.11). 13: Obtener la variable primal d. 14: Añadir la dirección extrema dl←d(con l=|DR|+ 1) al conjunto de direcciones extremas DR←DR∪ {l}. 15: Devolver la solución del problema original x=Pi∈ARλi·xi+ Pj∈DRµj·dj. A.3. Descomposición de Benders Consideremos el siguiente problema de programación lineal con variables complicantes, m´ın x∈Rn1,y∈Rn2c|·x+d|·y s.a: A·x+D·y≥a B·y≥b x≥0. (A.12) 116
APÉNDICE A. TÉCNICAS DE DESCOMPOSICIÓN Se considera que yes el vector de variables complicantes. Los vectores c∈Rn1yd∈Rn2son los coeficientes de las variables xeyen la función objetivo. Además las restricciones se dividen en dos conjuntos. En el primer conjunto de restricciones, A∈Rm1×n1yD∈Rm1×n2son las matrices de coeficientes de las variables xey, respectivamente, y a∈Rm1es el vector de términos independientes. El segundo conjunto de restricciones únicamente está definido a partir de las variables complicantes y, siendo B∈Rm2×n2su matriz de coeficientes y b∈Rm2su vector de términos independientes. El algoritmo de Benders [14] surge a principios de los años sesenta para la resolución de problemas con variables complicantes. Inicialmente, este algoritmo se crea para la resolución de problemas MILP, proponiendo la descomposición del problema teniendo en cuenta la naturaleza de sus variables. Con esta idea surgen las definiciones del problema maestro, construido con las variables enteras (complicantes), y del subproblema, construido con las variables continuas. De forma más general, las ideas del método de descomposición fueron rápidamente extendidas para otros tipos de problemas con variables complicantes. Esta flexibilidad del algoritmo de Benders es la principal razón de que haya sido aplicado exitosamente en multitud de problemas reales en distintos campos. Benders propuso un algoritmo iterativo que, como se ha dicho, se apoya en la resolución de dos problemas, el maestro y el subproblema. Estos no se resuelven de forma independiente. El problema maestro se comunica con el subproblema proponiéndole un valor de las variables complicantes mientras que el subproblema se comunica con el maestro a través de su solución dual, tal y como explicaremos más adelante tras una serie de preliminares necesarios. Para un valor fijo ¯yde yse define el subproblema como, m´ın x∈Rn1c|·x+d|·¯y s.a: A·x≥(a−D·¯y) x≥0, (A.13) 117
Jorge Rodríguez Veiga y su dual viene dado por, m´ax u∈Rm1(a−D·¯y)|·u+d|·¯y s.a: A|·u≤c u≥0. (A.14) Si definimos Y={y∈Rn2:B·y≥b}, entonces el problema original (A.12) se puede escribir empleando la formulación del subproblema (A.13) como m´ın y∈Ym´ın x∈Rn1{c|·x+d|·y:A·x+D·y≥a, x≥0}. Equivalentemente, por el Teorema de Dualidad Fuerte A.4 la definición se puede realizar a través de la formulación del subproblema dual (A.14), m´ın y∈Ym´ax u∈Rm1{(a−D·y)|·u+d|·y:A|·u≤c, u≥0}.(A.15) Consideremos ahora que el conjunto de puntos extremos de la región factible del subproblema dual es {u1, . . . , up}={ui∈Rm1:i∈A} y que el conjunto de sus direcciones extremas es {v1, . . . , vl}={vj∈ Rm1:j∈D}. Entonces, por el Teorema de Carathéodory A.1, se puede reformular el subproblema dual como, m´ax λ∈R|A|,µ∈R|D|X i∈A (a−D·¯y)|·ui·λi+X j∈D (a−D·¯y)|·vj·µj+d|·¯y s.a: X i∈A λi= 1 λi≥0,∀i∈A µj≥0,∀j∈D. El Corolario A.5 establece que para que el subproblema primal no sea infactible se ha de exigir que el subproblema dual sea acotado, i.e., (a−D·¯y)|·vj≤0para todo j∈D. Además, por el Teorema Fundamental de la Programación Lineal A.2, se sabe que existe un punto extremo de la región factible del subproblema dual que es óptimo. Por lo tanto, se 118
APÉNDICE A. TÉCNICAS DE DESCOMPOSICIÓN podrá reescribir como, m´ax i∈A(a−D·¯y)|·ui+d|·¯y (a−D·¯y)|·vj≤0,∀j∈D.(A.16) En consecuencia, podremos escribir el problema (A.15) como, m´ın y∈Ym´ax i∈A{(a−D·y)|·ui+d|·y: (a−D·y)|·vj≤0,∀j∈D}, o equivalentemente, empleando una variable auxiliar para modelar el máximo como, m´ın y∈Y, α∈Rα s.a: α≥d|·y+ (a−D·y)|·ui,∀i∈A (a−D·y)|·vj≤0,∀j∈D. (A.17) Al igual que el algoritmo de Dantzig-Wolfe, el algoritmo de Benders, propone relajar el problema (A.17) evitando considerar el conjunto total de puntos extremos y direcciones extremas de la región factible del subproblema dual. Partiendo de un subconjunto de los mismos, el algoritmo irá añadiendo nuevos puntos extremos o direcciones extremas. Sea por tanto AR⊂Aun subconjunto de puntos extremos y DR⊂Dun subconjunto de direcciones extremas. Entonces el problema maestro del algoritmo de Benders se define como, m´ın y∈Rn2,α∈Rα s.a: α≥d|·y+ (a−D·y)|·ui,∀i∈AR (a−D·y)|·vj≤0,∀j∈DR B·y≥b. (A.18) Las restricciones α≥d|·y+ (a−D·y)|·ui,∀i∈ARreciben el nombre de cortes de optimalidad, mientras que las restricciones (a−D· y)|·vj≤0,∀j∈DRson conocidas como cortes de factibilidad. El problema (A.18) puede tener un óptimo no finito aun cuando el problema original tenga un óptimo finito. Para evitar este comportamiento, se 119
Jorge Rodríguez Veiga pueden acotar las variables yyαconsiderando un valor Msuficientemente grande. De este modo, se podrían añadir al problema maestro las restricciones −M≤y≤My−M≤α. Por otra parte, si el subproblema dual es no acotado, se obtiene una dirección extrema resolviendo el siguiente problema, m´ax v∈Rm1(a−D·¯y)|·v s.a: A|·v≤0 ~ 1|·v≤1 v≥0, (A.19) donde ~ 1 = (1,...,1) ∈Rm1. La solución óptima del problema (A.19) es una dirección extrema de la región factible del subproblema dual. A continuación explicamos en líneas generales como actúa el algoritmo iterativo de Benders, explicando en particular la relación entre el problema maestro (A.18) y el subproblema dual (A.14). Se comienza resolviendo el problema maestro, obteniendo un valor ¯ypara la variable y así como una cota inferior para la función objetivo del problema original. Haciendo uso de este valor ¯yse resuelve el subproblema dual. Si el subproblema dual tiene una solución óptima finita, esta será un punto extremo de su región factible y el valor de la función objetivo será una cota superior para la función objetivo del problema original. Además se usa este punto extremo para incorporar al problema maestro un corte de optimalidad. Si el subproblema dual no es acotado, la resolución del problema (A.19) nos proporciona una nueva dirección extrema de la región factible. Dicha dirección se utiliza para incorporar en el problema maestro un corte de factibilidad. En cualquier caso, tras la incorporación del corte pertinente, se resuelve de nuevo el problema maestro y se repiten los pasos mientras que la cota inferior sea estrictamente menor que la superior. Nótese que el algoritmo finaliza en un número finito de iteraciones. Esto es así, pues en cada una de ellas se obtiene un punto o una dirección extrema de la región factible del subproblema dual (que es la misma en todas las iteraciones al no depender del valor ¯yconsiderado). Además, el número de puntos y direcciones extremas es finito. Por otra 120
APÉNDICE A. TÉCNICAS DE DESCOMPOSICIÓN parte, la incorporación de los cortes de factibilidad al problema maestro impide obtener dos veces la misma dirección extrema al resolver el problema (A.19). Finalmente, la obtención del mismo punto extremo en dos iteraciones implica el cumplimiento de la condición de parada del algoritmo. El Algoritmo 3 muestra el pseudocódigo del algoritmo de descomposición de Benders. Una descripción más detallada del mismo se puede encontrar en [35] y [45]. La flexibilidad del algoritmo de Benders permite que en muchas ocasiones sea habitual crear cortes específicos a partir del conocimiento del problema. Esto es algo habitual en problemas binarios o enteros. Ejemplos de ello son [50], [51] y [74]. Por último, indicar que existen diversas adaptaciones del algoritmo de Benders. Una de ellas consiste en adaptar el algoritmo para explotar la estructura del problema original permitiendo su descomposición en un problema maestro y varios subproblemas. Otro ejemplo consiste en su adaptación para la resolución de problemas no lineales [38]. Además, con el auge de la programación estocástica, es habitual encontrar en la literatura adaptaciones del algoritmo para la resolución de estos problemas. Estas adaptaciones son conocidas como L-Shape, para el caso en 2 etapas, o Nested Distance, para el caso multietapa [16]. 121
Jorge Rodríguez Veiga Algoritmo 3 Descomposición de Benders Input: Inicializar UB ← ∞,LB ← −∞, y los subconjuntos de puntos extremos y direcciones extremas AR=DR← ∅. 1: while LB < UB do 2: Resolver el problema maestro (A.18). 3: Obtener una solución óptima αey. 4: Establecer ¯y←yy actualizar LB ←α. 5: Resolver el subproblema dual (A.14). 6: if El subproblema es factible then 7: Obtener una solución primal xy dual uóptimas. 8: Añadir el punto extremo up←u(con p=|AR|+ 1) al conjunto de puntos extremos AR←AR∪ {p}y añadir el corte de optimalidad al problema maestro. 9: Actualizar UB ←c|·x+d|·¯y 10: else 11: Resolver el problema (A.19). 12: Obtener una solución óptima v(dirección extrema del subproblema dual). 13: Añadir la dirección extrema vl←v(con l=|DR|+ 1) al conjunto de direcciones extremas DR←DR∪ {l}y añadir el corte de factibilidad al problema maestro. 14: Devolver la solución (x, y). 122
Índice de figuras 1.1. Esquema de uso de la investigación operativa. . . . . . . . 5 3.1. Evolución histórica de los valores e índices de poder. . . . 23 4.1. Incendios con intervención de medios en 2018 en España. Copyright 2019 por Gobierno de España. Ministerio de Agricultura, Pesca y Alimentación [40]. Reimpreso con permiso. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 4.2. Ilustración del problema de Donovan y Rideout. . . . . . . 38 4.3. Ilustración del problema SARC. . . . . . . . . . . . . . . . . 39 4.4. Ilustración del problema AARV. . . . . . . . . . . . . . . . . 42 4.5. Ilustración del problema AAPR. . . . . . . . . . . . . . . . . 44 4.6. Diagrama de flujo para la contención eficiente de un incendio forestal. . . . . . . . . . . . . . . . . . . . . . . . . 47 123