Qual é o custo mínimo para conectar todas as ilhas?

Existe uma grelha de dimensão N x M . Algumas células são ilhas indicadas por '0' e as outras são água . Cada célula de água tem um número que denota o custo de uma ponte feita nessa célula. Você tem que encontrar o custo mínimo pelo qual todas as ilhas podem ser conectadas. Uma célula é conectada a outra célula se compartilhar uma aresta ou um vértice.

Qual algoritmo pode ser usado para resolver este problema? O que pode ser usado como uma abordagem de força bruta se os valores de N, M são muito pequenos, digamos NxM <= 100?

Exemplo : Na imagem fornecida, as células verdes indicam ilhas, as células azuis indicam água e as células azuis claras indicam as células nas quais uma ponte deve ser feita. Assim, para a imagem a seguir, a resposta será 17 .

Inicialmente pensei em marcar todas as ilhas como nós e conectar cada par de ilhas por uma ponte mais curta. Então, o problema poderia ser reduzido à árvore de abrangência mínima, mas nesta abordagem eu perdi o caso em que as arestas estão sobrepostas. Por exemplo , na imagem a seguir, a distância mais curta entre duas ilhas é 7 (marcada em amarelo), portanto, usando as Árvores de Extensão Mínima, a resposta seria 14 , mas a resposta deveria ser 11 (marcada em azul claro).

Para abordar este problema, eu usaria uma estrutura de programação inteira e definiria três conjuntos de variáveis ​​de decisão:

• x_ij : Uma variável indicadora binária para saber se construímos uma ponte no local da água (i, j).
• y_ijbcn : Um indicador binário para determinar se a localização da água (i, j) é a enésima localização ligando a ilha b à ilha c.
• l_bc : Uma variável indicadora binária para saber se as ilhas bec estão diretamente vinculadas (ou seja, você pode caminhar apenas em quadrados de pontes de b para c).

Para os custos de construção da ponte c_ij , o valor objetivo a minimizar é sum_ij c_ij * x_ij. Precisamos adicionar as seguintes restrições ao modelo:

• Precisamos garantir que as variáveis y_ijbcn sejam válidas. Sempre podemos chegar a um quadrado de água se construirmos uma ponte lá, portanto, y_ijbcn <= x_ijpara cada local de água (i, j). Além disso, y_ijbc1deve ser igual a 0 se (i, j) não faz fronteira com a ilha b. Finalmente, para n> 1, y_ijbcnsó pode ser usado se um local de água vizinho foi usado na etapa n-1. Definindo N(i, j)ser os quadrados de água vizinhos (i, j), isso é equivalente a y_ijbcn <= sum_{(l, m) in N(i, j)} y_lmbc(n-1).
• Precisamos garantir que as variáveis l_bc sejam definidas apenas se bec estiverem vinculados. Se definirmos I(c)como os locais que fazem fronteira com a ilha c, isso pode ser feito com l_bc <= sum_{(i, j) in I(c), n} y_ijbcn.
• Precisamos garantir que todas as ilhas estejam conectadas, direta ou indiretamente. Isso pode ser realizado da seguinte maneira: para cada subconjunto próprio não vazio S de ilhas, exigir que pelo menos uma ilha em S esteja ligada a pelo menos uma ilha no complemento de S, que chamaremos de S '. Em restrições, podemos implementar isso adicionando uma restrição para cada conjunto não vazio S de tamanho <= K / 2 (onde K é o número de ilhas) sum_{b in S} sum_{c in S'} l_bc >= 1,.

Para uma instância de problema com K ilhas, W quadrados de água e comprimento de caminho máximo especificado N, este é um modelo de programação inteira mista com O(K^2WN)variáveis ​​e O(K^2WN + 2^K)restrições. Obviamente, isso se tornará intratável à medida que o tamanho do problema se tornar grande, mas pode ser resolvido para os tamanhos de seu interesse. Para ter uma noção da escalabilidade, vou implementá-la em python usando o pacote pulp. Vamos primeiro começar com o mapa menor de 7 x 9 com 3 ilhas na parte inferior da questão:

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 ""

Isso leva 1,4 segundos para ser executado usando o solucionador padrão do pacote de celulose (o solucionador CBC) e produz a solução correta:

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

Em seguida, considere o problema completo no topo da questão, que é uma grade 13 x 14 com 7 ilhas:

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

Os solucionadores de MIP geralmente obtêm boas soluções com relativa rapidez e, então, gastam muito tempo tentando provar a otimização da solução. Usando o mesmo código de solucionador acima, o programa não é concluído em 30 minutos. No entanto, você pode fornecer um tempo limite para o solucionador para obter uma solução aproximada:

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

Isso resulta em uma solução com 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 melhorar a qualidade das soluções que você obtém, você pode usar um solucionador MIP comercial (gratuito se você estiver em uma instituição acadêmica e provavelmente não será gratuito caso você esteja). Por exemplo, aqui está o desempenho do Gurobi 6.0.4, novamente com um limite de tempo de 2 minutos (embora no log da solução tenhamos lido que o solucionador encontrou a melhor solução atual em 7 segundos):

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

Isso realmente encontra uma solução de valor objetivo 16, melhor do que o OP foi capaz de encontrar manualmente!

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 

Uma abordagem de força bruta, em 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

Em C ++, isso pode ser escrito 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);
}

Depois de fazer uma primeira chamada (presumo que você está transformando seus mapas 2d em arrays 1d para facilitar a cópia), bestCostconterá o custo da melhor resposta e bestconterá o padrão de pontes que o produz. No entanto, isso é extremamente lento.

Otimizações:

• Usando um "limite de ponte" e executando o algoritmo para aumentar o número máximo de pontes, você pode encontrar respostas mínimas sem explorar a árvore inteira. Encontrar uma resposta de 1 ponte, se existisse, seria O (nm) em vez de O (2 ^ nm) - uma melhoria drástica.
• Você pode evitar a pesquisa (interrompendo a recursão; isso também é chamado de "poda") depois de ter excedido bestCost, porque não faz sentido continuar procurando depois. Se não puder melhorar, não continue cavando.
• A poda acima funciona melhor se você olhar para os candidatos "bons" antes de olhar para os "ruins" (do jeito que está, as células são todas examinadas da esquerda para a direita, de cima para baixo). Uma boa heurística seria considerar as células que estão perto de vários componentes não conectados como de prioridade mais alta do que as células que não o são. No entanto, depois de adicionar heurísticas, sua pesquisa começa a se parecer com A * (e você também precisa de algum tipo de fila de prioridade).
• Pontes duplicadas e pontes para lugar nenhum devem ser evitadas. Qualquer ponte que não desconecte a rede da ilha, se removida, é redundante.

Um algoritmo de pesquisa geral como A * permite uma pesquisa muito mais rápida, embora encontrar melhores heurísticas não seja uma tarefa simples. Para uma abordagem mais específica do problema, usar os resultados existentes nas árvores Steiner , conforme sugerido por @Gassa, é o caminho a percorrer. Observe, entretanto, que o problema de construir árvores de Steiner em grades ortogonais é NP-Completo, de acordo com este artigo de Garey e Johnson .

Se "bom o suficiente" é o suficiente, um algoritmo genético provavelmente pode encontrar soluções aceitáveis ​​rapidamente, desde que você adicione algumas heurísticas importantes quanto ao posicionamento preferido da ponte.

Este problema é uma variante da árvore de Steiner chamada árvore de Steiner ponderada por nós , especializada em uma determinada classe de grafos. De forma compacta, a árvore de Steiner com peso de nó é, dado um grafo não direcionado com peso de nó onde alguns nós são terminais, encontrar o conjunto mais barato de nós, incluindo todos os terminais que induzem um subgrafo conectado. Infelizmente, não consigo encontrar nenhum solucionador em algumas pesquisas superficiais.

Para formular como um programa inteiro, faça uma variável 0-1 para cada nó não terminal, então para todos os subconjuntos de nós não terminais cuja remoção do gráfico inicial desconecta dois terminais, requer que a soma das variáveis ​​no subconjunto esteja em pelo menos 1. Isso induz muitas restrições, então você terá que aplicá-las preguiçosamente, usando um algoritmo eficiente para conectividade de nó (fluxo máximo, basicamente) para detectar uma restrição violada ao máximo. Desculpem a falta de detalhes, mas vai ser difícil implementar, mesmo se você já estiver familiarizado com a programação de inteiros.

Dado que esse problema ocorre em uma grade e você tem parâmetros bem definidos, eu abordaria o problema com a eliminação sistemática do espaço do problema por meio da criação de uma árvore geradora mínima. Fazendo isso, faz sentido para mim abordar esse problema com o algoritmo de Prim.

Infelizmente, agora você se depara com o problema de abstrair a grade para criar um conjunto de nós e arestas ... logo, o verdadeiro problema deste post é como converter minha grade nxm em {V} e {E}?

Este processo de conversão é, à primeira vista, provavelmente NP-Difícil devido ao grande número de combinações possíveis (suponha que todos os custos hidroviários sejam idênticos). Para lidar com instâncias onde os caminhos se sobrepõem, você deve considerar fazer uma ilha virtual.

Quando isso for feito, execute o Algoritmo de Prim e você deve chegar à solução ótima.

Não acredito que a Programação Dinâmica possa ser executada de forma eficaz aqui porque não existe um princípio observável de otimização. Se encontrarmos o custo mínimo entre duas ilhas, isso não significa necessariamente que podemos encontrar o custo mínimo entre essas duas e a terceira ilha, ou outro subconjunto de ilhas que será (por minha definição encontrar o MST via Prim) conectado.

Se você deseja que o código (pseudo ou outro) converta sua grade em um conjunto de {V} e {E}, envie-me uma mensagem privada e irei ver como unir uma implementação.