Encontrando o par de pontos mais próximo

Declaração do problema

Dados $n$ pontos no plano. Cada ponto $p_i$ é definido por suas coordenadas $(x_i,y_i)$. É necessário encontrar entre eles dois desses pontos, de modo que a distância entre eles seja mínima:

$$ \min_{\scriptstyle i, j=0 \ldots n-1,\atop \scriptstyle i \neq j } \rho (p_i, p_j). $$

Tomamos as distâncias euclidianas usuais:

$$ \rho (p_i,p_j) = \sqrt{(x_i-x_j)^2 + (y_i-y_j)^2} .$$

O algoritmo trivial - iterando sobre todos os pares e calculando a distância para cada um - funciona em $O(n^2)$.

O algoritmo executado em tempo $O(n \log n)$ é descrito abaixo. Esse algoritmo foi proposto por Preparata em 1975. Preparata e Shamos também mostraram que esse algoritmo é ideal no modelo de uma árvore de decisões.

Algoritmo

Construímos um algoritmo de acordo com o esquema geral dos algoritmos de dividir e conquistar: o algoritmo é projetado como uma função recursiva, para a qual passamos um conjunto de pontos; essa função recursiva divide esse conjunto pela metade, chama-se recursivamente em cada metade e executa algumas operações para combinar as respostas. A operação de combinação consiste em detectar os casos em que um ponto da solução ideal caiu na metade e o outro no outro (nesse caso, chamadas recursivas de cada uma das metades não podem detectar esse par separadamente). A principal dificuldade, como sempre no caso de algoritmos de dividir e conquistar, reside na implementação efetiva do estágio de mesclagem. Se um conjunto de $n$ pontos for passado para a função recursiva, o estágio de mesclagem não deve funcionar mais do que O(n), então o comportamento assintótico de todo o algoritmo $T(n)$ será encontrado na equação:

$$T(n) = 2T(n/2) + O(n).$$

A solução para esta equação é $T(n) = O(n \log n).$

Então, prosseguimos para a construção do algoritmo. Para chegar a uma implementação eficaz do estágio de mesclagem no futuro, dividiremos o conjunto de pontos em dois subconjuntos, de acordo com suas coordenadas $x$: Na verdade, "desenhamos" uma linha vertical dividindo o conjunto de pontos em dois subconjuntos aproximadamente do mesmo tamanho. É conveniente fazer essa partição da seguinte maneira: ordenamos os pontos da maneira padrão como pares de números, ou seja:

$$p_i < p_j \Longleftrightarrow (x_i < x_j) \lor \Big(\left(x_i = x_j\right) \wedge \left(y_i < y_j \right) \Big) $$

Em seguida, pegamos o ponto do meio depois de ordenar $p_m (m = \lfloor n/2 \rfloor)$, e todos os pontos antes dele e o próprio $p_m$ são atribuídos à primeira metade e todos os pontos depois dele - segunda metade:

$$A_1 = {p_i \ | \ i = 0 \ldots m }$$ $$A_2 = {p_i \ | \ i = m + 1 \ldots n-1 }.$$

Agora, chamando recursivamente cada um dos conjuntos $A_1$ e $A_2$, encontraremos as respostas $h_1$ e $h_2$ para cada uma das metades. E então pegamos o melhor deles $h = \min(h_1, h_2)$.

Agora precisamos fazer um estágio de mesclagem , ou seja, tentamos encontrar esses pares de pontos, para os quais a distância entre eles seja menor que $h$ e um ponto estar em $A_1$ e o outro em $A_2$. É suficiente considerar apenas os pontos que são separados da linha vertical por uma distância inferior a $h$, ou seja, o conjunto $B$ dos pontos considerados nesta etapa é igual a:

$$B = { p_i\ | \ | x_i - x_m\ | < h }.$$

Para cada ponto no conjunto $B$, tentamos encontrar os pontos que estão mais próximos do que $h$. Por exemplo, é suficiente considerar apenas os pontos cuja coordenada $y$ difere em não mais que $h$. Além disso, não faz sentido considerar os pontos cuja coordenada $y$ é maior que a coordenada $y$ do ponto atual. Assim, para cada ponto $p_i$ definimos o conjunto de pontos considerados $C(p_i)$ da seguinte maneira:

$$C(p_i) = { p_j\ |\ p_j \in B,\ \ y_i - h < y_j \le y_i }.$$

Se ordenarmos os pontos do conjunto $B$ pela coordenada $y$, será mais fácil encontrar $C(p_i)$: esses são vários pontos na linha seguidos do ponto $p_i$.

Portanto, na nova notação, o estágio de mesclagem se parece com isso: construa um conjunto $B$, ordene os pontos nele pela coordenada $y$, e então para cada ponto $p_i \in B$ considere todos os pontos $p_j \in C(p_i)$, e para cada par $(p_i,p_j)$ calcule a distância e compare com a melhor distância atual.

À primeira vista, este ainda é um algoritmo não ideal: parece que os tamanhos dos conjuntos $C(p_i)$ serão da ordem de $n$, e o comportamento assintótico necessário não funcionará. No entanto, surpreendentemente, pode-se provar que o tamanho de cada um dos conjuntos $C(p_i)$ é uma quantidade $O(1)$, ou seja, não excede alguma pequena constante, independentemente dos próprios pontos. A prova desse fato é dada na próxima seção.

Finalmente, prestamos atenção à ordenação(sorting) que o algoritmo acima contém: primeiro, ordenamos pelos pares $(x, y)$, e, em seguida, ordenamos os elementos do conjunto $B$ por $y$. De fato, esses duas ordenações dentro da função recursiva podem ser eliminadas (caso contrário, não atingiríamos a estimativa $O(n)$ para o estágio de mesclagem, e o comportamento assintótico do algoritmo iria ser $O(n \log^2 n)$). É fácil se livrar da primeira ordenação - basta fazer essa ordenação antes de iniciar a recursão: afinal, os próprios elementos não mudam dentro da recursão, portanto, não há necessidade de ordenar novamente. Com a segunda ordenação é um pouco mais difícil, executá-la anteriormente não funcionará. Mas, lembrando o tipo de mesclagem, que também funciona com o princípio de dividir e conquistar, podemos simplesmente incorporar esse tipo em nossa recursão. Deixe a recursão, pegando algum conjunto de pontos (como lembramos, ordenados pelos pares $(x, y)$), retornando o mesmo conjunto, porém, ordenados pela coordenada $y$. Para fazer isso, simplesmente basta mesclar (em $O(n)$) os dois resultados retornados por chamadas recursivas. Isso resultará em um conjunto ordenado pela coordenada $y$.

Análise assintótica

Para mostrar que o algoritmo acima é realmente executado em $O(n \log n)$, precisamos provar o seguinte fato: $|C(p_i)| = O(1)$.

Então, vamos considerar algum ponto $p_i$; lembre-se de que o conjunto $C(p_i)$ é um conjunto de pontos cuja coordenada $y$ está no segmento $[y_i-h; y_i]$, e, além disso, ao longo da coordenada $x$, do próprio ponto $p_i$ e de todos os pontos do conjunto $C(p_i)$ estão dentro do limite(largura) $2h$. Em outras palavras, os pontos que estamos considerando $p_i$ e $C(p_i)$ estão em um retângulo de tamanho $2h \times h$.

Nossa tarefa é estimar o número máximo de pontos que podem estar nesse retângulo $2h \times h$; portanto, estimamos o tamanho máximo do conjunto $C(p_i)$. Ao mesmo tempo, ao avaliar, não devemos esquecer que pode haver pontos repetidos.

Lembre-se de que $h$ foi obtido a partir dos resultados de duas chamadas recursivas - nos conjuntos $A_1$ e $A_2$, e $A_1$ contém pontos à esquerda da linha de partição e parcialmente nela, $A_2$ contém os pontos restantes da linha de partição e pontos para a direita dela. Para qualquer par de pontos de $A_1$, também de $A_2$, a distância não pode ser menos que $h$ — caso contrário, isso significaria uma operação incorreta da função recursiva.

Para estimar o número máximo de pontos no retângulo $2h \times h$ dividimos em dois quadrados $h \times h$, o primeiro quadrado inclui todos os pontos $C(p_i) \cap A_1$, e o segundo contém todos os outros, ou seja, $C(p_i) \cap A_2$. Resulta das considerações acima que em cada um desses quadrados a distância entre dois pontos é de pelo menos $h$.

Mostramos que existem no máximo quatro pontos em cada quadrado. Por exemplo, isso pode ser feito da seguinte maneira: divida o quadrado em $4$ sub-quadrados com lados $h/2$. Então não pode haver mais de um ponto em cada um desses sub-quadrados (já que até a diagonal é igual a $h / \sqrt{2}$, o que é menor que $h$). Portanto, não pode haver mais que $4$ pontos em todo o quadrado.

Portanto, provamos que em um retângulo $2h \times h$ não pode haver mais que $4 \cdot 2 = 8$ pontos, e, portanto, o tamanho do conjunto $C(p_i)$ não pode exceder $7$, como requerido.

Implementação

Introduzimos uma estrutura de dados para armazenar um ponto (suas coordenadas e um número) e operadores de comparação necessários para dois tipos de ordenação:

struct pt {
    int x, y, id;
};

struct cmp_x {
    bool operator()(const pt & a, const pt & b) const {
        return a.x < b.x || (a.x == b.x && a.y < b.y);
    }
};

struct cmp_y {
    bool operator()(const pt & a, const pt & b) const {
        return a.y < b.y;
    }
};

int n;
vector<pt> a;

Para uma implementação conveniente da recursão, apresentamos uma função auxiliar upd_ans(), que calcula a distância entre dois pontos e verifica se é melhor que a resposta atual:

double mindist;
pair<int, int> best_pair;

void upd_ans(const pt & a, const pt & b) {
    double dist = sqrt((a.x - b.x)*(a.x - b.x) + (a.y - b.y)*(a.y - b.y));
    if (dist < mindist) {
        mindist = dist;
        best_pair = {a.id, b.id};
    }
}

Finalmente, a implementação da recursão em si. Supõe-se que, antes de chamá-la, a array $a[]$ ja foi ordenada pela coordenada $x$. Na recursão, passamos apenas dois ponteiros $l, r$, o que indica que ele deve procurar a resposta para $a[l \ldots r)$. Se a distância entre $r$ e $l$ for muito pequena, a recursão deve ser interrompida e deve ser executado um algoritmo trivial para encontrar o par mais próximo e, em seguida, ordenar a subarray pela coordenada $y$.

Para mesclar dois conjuntos de pontos, recebidos de chamadas recursivas, em um (ordenados por $y$), usamos a função STL $merge()$, e criamos um buffer auxiliar $t[]$(um para todas as chamadas recursivas). (Usar inplace_merge () é impraticável porque geralmente não funciona em tempo linear)

Finalmente, o conjunto $B$ é armazenado na mesma array $t$.

vector<pt> t;

void rec(int l, int r) {
    if (r - l <= 3) {
        for (int i = l; i < r; ++i) {
            for (int j = i + 1; j < r; ++j) {
                upd_ans(a[i], a[j]);
            }
        }
        sort(a.begin() + l, a.begin() + r, cmp_y());
        return;
    }

    int m = (l + r) >> 1;
    int midx = a[m].x;
    rec(l, m);
    rec(m, r);

    merge(a.begin() + l, a.begin() + m, a.begin() + m, a.begin() + r, t.begin(), cmp_y());
    copy(t.begin(), t.begin() + r - l, a.begin() + l);

    int tsz = 0;
    for (int i = l; i < r; ++i) {
        if (abs(a[i].x - midx) < mindist) {
            for (int j = tsz - 1; j >= 0 && a[i].y - t[j].y < mindist; --j)
                upd_ans(a[i], t[j]);
            t[tsz++] = a[i];
        }
    }
}

A propósito, se todas as coordenadas forem inteiras, no momento da recursão, você não poderá mover para valores fracionários e armazenar em $mindist$ o quadrado da distância mínima.

No programa principal("main()"), a recursão deve ser chamada da seguinte maneira:

t.resize(n);
sort(a.begin(), a.end(), cmp_x());
mindist = 1E20;
rec(0, n);

Generalização: encontrando um triângulo com perímetro mínimo

O algoritmo descrito acima é curiosamente generalizado para esse problema: entre um determinado conjunto de pontos, escolha três pontos diferentes para que a soma das distâncias de pares entre eles seja a menor.

De fato, para resolver esse problema, o algoritmo permanece o mesmo: dividimos o campo em duas metades da linha vertical, chamamos a solução recursivamente nas duas metades, escolhemos o mínimo $minper$ a partir dos perímetros encontrados, construimos uma "faixa" com o "espessura" de $minper / 2$, e iteramos por todos os triângulos que podem melhorar a resposta. (Observe que o triângulo com perímetro $\le minper$ tem o maior tamanho $\le minper / 2$.)

Problemas