algorithm dynamic-programming mathematical-optimization linear-programming heuristics

algorithm - ¿Cuál es el costo mínimo para conectar todas las islas?



dynamic-programming mathematical-optimization (4)

Dado que este problema tiene lugar en una cuadrícula y tiene parámetros bien definidos, abordaría el problema con la eliminación sistemática del espacio del problema mediante la creación de un árbol de expansión mínimo. Al hacerlo, tiene sentido para mí si abordas este problema con el algoritmo de Prim.

Desafortunadamente, ahora se encuentra con el problema de abstraer la cuadrícula para crear un conjunto de nodos y bordes ... ergo, el verdadero problema de esta publicación es cómo convierto mi cuadrícula nxm en {V} y {E}.

Este proceso de conversión es, de un vistazo, probablemente NP-Hard debido a la gran cantidad de combinaciones posibles (suponga que todos los costos de las vías navegables son idénticos). Para manejar instancias donde las rutas se superponen, debe considerar hacer una isla virtual.

Cuando haya terminado, ejecute el Algoritmo de Prim y debería llegar a la solución óptima.

No creo que la programación dinámica se pueda ejecutar de manera efectiva aquí porque no hay un principio observable de optimización. Si encontramos el costo mínimo entre dos islas, eso no significa necesariamente que podamos encontrar el costo mínimo entre esas dos y la tercera isla, u otro subconjunto de islas que será (según mi definición, encontrar el MST a través de Prim) conectado.

Si desea que el código (pseudo o de otro tipo) convierta su grilla en un conjunto de {V} y {E}, envíeme un mensaje privado y buscaré unir una implementación.

Hay una cuadrícula de tamaño N x M. Algunas celdas son islas denotadas por ''0'' y las otras son agua . Cada celda de agua tiene un número que indica el costo de un puente hecho en esa celda. Debe encontrar el costo mínimo por el cual se pueden conectar todas las islas. Una celda está conectada a otra celda si comparte un borde o un vértice.

¿Qué algoritmo se puede usar para resolver este problema?
Editar: ¿Qué se puede usar como un enfoque de fuerza bruta si los valores de N, M son muy pequeños, digamos NxM <= 100?

Ejemplo : en la imagen dada, las celdas verdes indican islas, las celdas azules indican agua y las celdas azules claras indican las celdas sobre las que se debe hacer un puente. Por lo tanto, para la siguiente imagen, la respuesta será 17 .

Inicialmente pensé en marcar todas las islas como nodos y conectar cada par de islas por un puente más corto. Entonces, el problema podría reducirse a un árbol de expansión mínimo, pero en este enfoque omití el caso en el que los bordes se superponen. Por ejemplo , en la siguiente imagen, la distancia más corta entre dos islas es 7 (marcada en amarillo), por lo que al usar Árboles de expansión mínima la respuesta sería 14 , pero la respuesta debería ser 11 (marcada en azul claro).

Este problema es una variante del árbol de Steiner llamado árbol de Steiner ponderado por nodo , especializado en una cierta clase de gráficos. Compactamente, el árbol Steiner ponderado por nodo es, dado un gráfico no dirigido ponderado por nodo donde algunos nodos son terminales, encuentre el conjunto de nodos más barato, incluidos todos los terminales que inducen un subgrafo conectado. Lamentablemente, parece que no puedo encontrar ningún solucionador en algunas búsquedas superficiales.

Para formular como un programa entero, haga una variable 0-1 para cada nodo no terminal, luego para todos los subconjuntos de nodos terminales cuya eliminación del gráfico inicial desconecte dos terminales, requiera que la suma de las variables en el subconjunto esté en menos 1. Esto induce demasiadas restricciones, por lo que tendrá que aplicarlas perezosamente, utilizando un algoritmo eficiente para la conectividad de nodos (flujo máximo, básicamente) para detectar una restricción máximamente violada. Perdón por la falta de detalles, pero será difícil de implementar, incluso si ya está familiarizado con la programación de enteros.

Para abordar este problema, usaría un marco de programación entero y definiría tres conjuntos de variables de decisión:

  • x_ij : una variable indicadora binaria para determinar si construimos un puente en la ubicación del agua (i, j).
  • y_ijbcn : un indicador binario para determinar si la ubicación del agua (i, j) es la enésima ubicación que une la isla b con la isla c.
  • l_bc : una variable indicadora binaria para determinar si las islas byc están directamente vinculadas (también conocido como solo se puede caminar en los cuadrados del puente de b a c).

Para los costos de construcción de puentes c_ij , el valor objetivo para minimizar es sum_ij c_ij * x_ij . Necesitamos agregar las siguientes restricciones al modelo:

  • Necesitamos asegurarnos de que las variables y_ijbcn sean válidas. Siempre podemos alcanzar un cuadrado de agua si construimos un puente allí, entonces y_ijbcn <= x_ij para cada ubicación de agua (i, j). Además, y_ijbc1 debe ser igual a 0 si (i, j) no limita con la isla b. Finalmente, para n> 1, y_ijbcn solo se puede usar si se usó una ubicación de agua vecina en el paso n-1. Definiendo N(i, j) como los cuadrados de agua vecinos (i, j), esto es equivalente a y_ijbcn <= sum_{(l, m) in N(i, j)} y_lmbc(n-1) .
  • Necesitamos asegurarnos de que las variables l_bc solo estén establecidas si byc están vinculadas. Si definimos I(c) como las ubicaciones que bordean la isla c, esto se puede lograr con l_bc <= sum_{(i, j) in I(c), n} y_ijbcn .
  • Necesitamos asegurarnos de que todas las islas estén conectadas, ya sea directa o indirectamente. Esto se puede lograr de la siguiente manera: para cada subconjunto de islas S no vacío, se requiere que al menos una isla en S esté vinculada a al menos una isla en el complemento de S, que llamaremos S ''. En restricciones, podemos implementar esto agregando una restricción para cada conjunto S no vacío de tamaño <= K / 2 (donde K es el número de islas), sum_{b in S} sum_{c in S''} l_bc >= 1 .

Para una instancia problemática con islas K, cuadrados de agua W y longitud máxima de ruta especificada N, este es un modelo de programación de enteros mixtos con variables O(K^2WN + 2^K) restricciones O(K^2WN + 2^K) . Obviamente, esto se volverá intratable a medida que el tamaño del problema se agrande, pero puede resolverse para los tamaños que le interesan. Para tener una idea de la escalabilidad, la implementaré en python usando el paquete pulp. Primero comencemos con el mapa más pequeño de 7 x 9 con 3 islas en la parte inferior de la pregunta:

import itertools import pulp water = {(0, 2): 2.0, (0, 3): 1.0, (0, 4): 1.0, (0, 5): 1.0, (0, 6): 2.0, (1, 0): 2.0, (1, 1): 9.0, (1, 2): 1.0, (1, 3): 9.0, (1, 4): 9.0, (1, 5): 9.0, (1, 6): 1.0, (1, 7): 9.0, (1, 8): 2.0, (2, 0): 1.0, (2, 1): 9.0, (2, 2): 9.0, (2, 3): 1.0, (2, 4): 9.0, (2, 5): 1.0, (2, 6): 9.0, (2, 7): 9.0, (2, 8): 1.0, (3, 0): 9.0, (3, 1): 1.0, (3, 2): 9.0, (3, 3): 9.0, (3, 4): 5.0, (3, 5): 9.0, (3, 6): 9.0, (3, 7): 1.0, (3, 8): 9.0, (4, 0): 9.0, (4, 1): 9.0, (4, 2): 1.0, (4, 3): 9.0, (4, 4): 1.0, (4, 5): 9.0, (4, 6): 1.0, (4, 7): 9.0, (4, 8): 9.0, (5, 0): 9.0, (5, 1): 9.0, (5, 2): 9.0, (5, 3): 2.0, (5, 4): 1.0, (5, 5): 2.0, (5, 6): 9.0, (5, 7): 9.0, (5, 8): 9.0, (6, 0): 9.0, (6, 1): 9.0, (6, 2): 9.0, (6, 6): 9.0, (6, 7): 9.0, (6, 8): 9.0} islands = {0: [(0, 0), (0, 1)], 1: [(0, 7), (0, 8)], 2: [(6, 3), (6, 4), (6, 5)]} N = 6 # Island borders iborders = {} for k in islands: iborders[k] = {} for i, j in islands[k]: for dx in [-1, 0, 1]: for dy in [-1, 0, 1]: if (i+dx, j+dy) in water: iborders[k][(i+dx, j+dy)] = True # Create models with specified variables x = pulp.LpVariable.dicts("x", water.keys(), lowBound=0, upBound=1, cat=pulp.LpInteger) pairs = [(b, c) for b in islands for c in islands if b < c] yvals = [] for i, j in water: for b, c in pairs: for n in range(N): yvals.append((i, j, b, c, n)) y = pulp.LpVariable.dicts("y", yvals, lowBound=0, upBound=1) l = pulp.LpVariable.dicts("l", pairs, lowBound=0, upBound=1) mod = pulp.LpProblem("Islands", pulp.LpMinimize) # Objective mod += sum([water[k] * x[k] for k in water]) # Valid y for k in yvals: i, j, b, c, n = k mod += y[k] <= x[(i, j)] if n == 0 and not (i, j) in iborders[b]: mod += y[k] == 0 elif n > 0: mod += y[k] <= sum([y[(i+dx, j+dy, b, c, n-1)] for dx in [-1, 0, 1] for dy in [-1, 0, 1] if (i+dx, j+dy) in water]) # Valid l for b, c in pairs: mod += l[(b, c)] <= sum([y[(i, j, B, C, n)] for i, j, B, C, n in yvals if (i, j) in iborders[c] and B==b and C==c]) # All islands connected (directly or indirectly) ikeys = islands.keys() for size in range(1, len(ikeys)/2+1): for S in itertools.combinations(ikeys, size): thisSubset = {m: True for m in S} Sprime = [m for m in ikeys if not m in thisSubset] mod += sum([l[(min(b, c), max(b, c))] for b in S for c in Sprime]) >= 1 # Solve and output mod.solve() for row in range(min([m[0] for m in water]), max([m[0] for m in water])+1): for col in range(min([m[1] for m in water]), max([m[1] for m in water])+1): if (row, col) in water: if x[(row, col)].value() > 0.999: print "B", else: print "-", else: print "I", print ""

Esto tarda 1,4 segundos en ejecutarse utilizando el solucionador predeterminado del paquete de pulpa (el solucionador CBC) y genera la solución correcta:

I I - - - - - I I - - B - - - B - - - - - B - B - - - - - - - B - - - - - - - - B - - - - - - - - B - - - - - - - I I I - - -

A continuación, considere el problema completo en la parte superior de la pregunta, que es una cuadrícula de 13 x 14 con 7 islas:

water = {(i, j): 1.0 for i in range(13) for j in range(14)} islands = {0: [(0, 0), (0, 1), (1, 0), (1, 1), (2, 0), (2, 1)], 1: [(9, 0), (9, 1), (10, 0), (10, 1), (10, 2), (11, 0), (11, 1), (11, 2), (12, 0)], 2: [(0, 7), (0, 8), (1, 7), (1, 8), (2, 7)], 3: [(7, 7), (8, 6), (8, 7), (8, 8), (9, 7)], 4: [(0, 11), (0, 12), (0, 13), (1, 12)], 5: [(4, 10), (4, 11), (5, 10), (5, 11)], 6: [(11, 8), (11, 9), (11, 13), (12, 8), (12, 9), (12, 10), (12, 11), (12, 12), (12, 13)]} for k in islands: for i, j in islands[k]: del water[(i, j)] for i, j in [(10, 7), (10, 8), (10, 9), (10, 10), (10, 11), (10, 12), (11, 7), (12, 7)]: water[(i, j)] = 20.0 N = 7

Los solucionadores de MIP a menudo obtienen buenas soluciones con relativa rapidez y luego pasan mucho tiempo tratando de demostrar la optimización de la solución. Usando el mismo código de solución que el anterior, el programa no se completa en 30 minutos. Sin embargo, puede proporcionar un tiempo de espera al solucionador para obtener una solución aproximada:

mod.solve(pulp.solvers.PULP_CBC_CMD(maxSeconds=120))

Esto produce una solución con valor objetivo 17:

I I - - - - - I I - - I I I I I - - - - - I I - - - I - I I - - - - - I - B - B - - - - B - - - B - - - B - - - - - - B - B - - - - I I - - - - - - B - - - - - I I - - - - - - - B - - - - - B - - - - - - - B - I - - - - B - - - - - B - I I I - - B - - I I - B - - - I - - - - B - I I I - - - - - - - - - - B I I I - - - - - I I - - - I I - - - - - - - I I I I I I

Para mejorar la calidad de las soluciones que obtiene, puede utilizar un solucionador de MIP comercial (esto es gratuito si se encuentra en una institución académica y probablemente no lo sea de otra manera). Por ejemplo, aquí está el rendimiento de Gurobi 6.0.4, nuevamente con un límite de tiempo de 2 minutos (aunque del registro de soluciones leemos que el solucionador encontró la mejor solución actual en 7 segundos):

mod.solve(pulp.solvers.GUROBI(timeLimit=120))

¡Esto realmente encuentra una solución de valor objetivo 16, una mejor que la OP pudo encontrar a mano!

I I - - - - - I I - - I I I I I - - - - - I I - - - I - I I - - - - - I - B - B - - - - B - - - - - - - B - - - - - - B - - - - - - I I - - - - - - B - - - - - I I - - - - - - - B - - B B - - - - - - - - - B - I - - B - - - - - - - B - I I I - - B - - I I - B - - - I - - - - B - I I I - - - - - - - - - - B I I I - - - - - I I - - - I I - - - - - - - I I I I I I

Un enfoque de fuerza bruta, en pseudocódigo:

start with a horrible "best" answer given an nxm map, try all 2^(n*m) combinations of bridge/no-bridge for each cell if the result is connected, and better than previous best, store it return best

En C ++, esto podría escribirse como

// map = linearized map; map[x*n + y] is the equivalent of map2d[y][x] // nm = n*m // bridged = true if bridge there, false if not. Also linearized // nBridged = depth of recursion (= current bridge being considered) // cost = total cost of bridges in ''bridged'' // best, bestCost = best answer so far. Initialized to "horrible" void findBestBridges(char map[], int nm, bool bridged[], int nBridged, int cost, bool best[], int &bestCost) { if (nBridged == nm) { if (connected(map, nm, bridged) && cost < bestCost) { memcpy(best, bridged, nBridged); bestCost = best; } return; } if (map[nBridged] != 0) { // try with a bridge there bridged[nBridged] = true; cost += map[nBridged]; // see how it turns out findBestBridges(map, nm, bridged, nBridged+1, cost, best, bestCost); // remove bridge for further recursion bridged[nBridged] = false; cost -= map[nBridged]; } // and try without a bridge there findBestBridges(map, nm, bridged, nBridged+1, cost, best, bestCost); }

Después de hacer una primera llamada (supongo que está transformando sus mapas 2D en matrices 1d para facilitar la copia), bestCost contendrá el costo de la mejor respuesta y el best contendrá el patrón de puentes que lo genera. Esto es, sin embargo, extremadamente lento.

Optimizaciones:

  • Al usar un "límite de puente" y ejecutar el algoritmo para aumentar el número máximo de puentes, puede encontrar respuestas mínimas sin explorar todo el árbol. Encontrar una respuesta de 1 puente, si existiera, sería O (nm) en lugar de O (2 ^ nm), una mejora drástica.
  • Puede evitar la búsqueda (deteniendo la recursividad; esto también se llama "poda") una vez que haya excedido bestCost , porque no tiene sentido seguir buscando después. Si no puede mejorar, no siga cavando.
  • La poda anterior funciona mejor si observa a los candidatos "buenos" antes de mirar a los "malos" (tal como están, las celdas se analizan en orden de izquierda a derecha, de arriba a abajo). Una buena heurística sería considerar las celdas que están cerca de varios componentes no conectados como de mayor prioridad que las celdas que no lo están. Sin embargo, una vez que agrega heurística, su búsqueda comienza a parecerse a A* (y también necesita algún tipo de cola prioritaria).
  • Los puentes duplicados y los puentes a ninguna parte deben ser evitados. Cualquier puente que no desconecte la red de la isla si se elimina es redundante.

Un algoritmo de búsqueda general como A* permite una búsqueda mucho más rápida, aunque encontrar mejores heurísticas no es una tarea sencilla. Para un enfoque más específico del problema, utilizar los resultados existentes en los árboles Steiner , como lo sugiere @Gassa, es el camino a seguir. Sin embargo, tenga en cuenta que el problema de construir árboles Steiner en cuadrículas ortogonales es NP-Complete, según este documento de Garey y Johnson .

Si "lo suficientemente bueno" es suficiente, un algoritmo genético probablemente pueda encontrar soluciones aceptables rápidamente, siempre que agregue algunas heurísticas clave en cuanto a la colocación preferida del puente.