#️⃣ Haszowanie przestrzenne: wyszukiwanie najbliższych sąsiadów O(1) dla symulacji cząstek

Naiwny problem O(n²)

Dla n cząstek naiwne znalezienie wszystkich par w odległości h wymaga n(n-1)/2 sprawdzeń odległości. Dla n = 1000 to 499 500 sprawdzeń. Dla n = 10 000 to 50 milionów. Przy 60 fps każda klatka ma około 16 ms — 50 milionów obliczeń odległości (każde to pierwiastek kwadratowy lub co najmniej trzy mnożenia i dodawanie) zajmie setki milisekund. Symulacja zawiesza się poniżej 10 fps.

Kluczową obserwacją jest to, że cząstki oddziałują tylko w promieniu h. Jeśli potrafimy szybko znaleźć, które cząstki znajdują się w odległości h od danej cząstki bez sprawdzania wszystkich pozostałych, możemy zredukować pracę z O(n²) do O(n·k), gdzie k to średnia liczba sąsiadów — zazwyczaj mała stała (5–50 dla typowych symulacji SPH).

Podejście siatkowe: podział przestrzenny

Podziel obszar symulacji na siatkę o rozmiarze komórki równym promieniowi oddziaływania h. Cząstka w pozycji (x, y) należy do komórki (⌊x/h⌋, ⌊y/h⌋). Ponieważ cząstki oddziałują tylko w odległości h, każda cząstka oddziałująca z cząstką w komórce (cx, cy) musi znajdować się w jednej z 9 sąsiednich komórek (w 2D) — samej komórce plus 8 otaczających. W 3D jest to 27 komórek.

Jeśli siatka ma G = (W/h) × (H/h) komórek (dla pudełka symulacji W × H), a n cząstek jest rozmieszczonych mniej więcej równomiernie, każda komórka zawiera średnio około n/G cząstek. Przeszukanie 9 komórek zawierających po n/G cząstek kosztuje O(9 · n/G) = O(n · h²/(WH)) pracy — proporcjonalnie do ułamka powierzchni pokrytej przez promień oddziaływania.

Haszowanie przestrzenne: bez alokowania siatki

Siatka o stałym rozmiarze wymaga zaalokowania W/h × H/h komórek — dobrze dla znanych rozmiarów domeny, ale marnotrawne, jeśli cząstki skupiają się w małym regionie lub domena jest duża. Haszowanie przestrzenne unika jawnej siatki, haszując współrzędne komórek do tablicy o stałym rozmiarze:

class SpatialHash {
  constructor(cellSize, tableSize) {
    this.cellSize = cellSize;
    this.tableSize = tableSize;
    this.table = new Array(tableSize).fill(null).map(() => []);
  }

  hash(cx, cy) {
    // Duże liczby pierwsze, aby równomiernie rozłożyć współrzędne siatki
    return ((cx * 92837111) ^ (cy * 689287499)) % this.tableSize;
  }

  cellCoords(x, y) {
    return [Math.floor(x / this.cellSize), Math.floor(y / this.cellSize)];
  }

  insert(particle) {
    const [cx, cy] = this.cellCoords(particle.x, particle.y);
    const key = this.hash(cx, cy);
    this.table[key].push(particle);
  }

  query(x, y) {
    const [cx, cy] = this.cellCoords(x, y);
    const neighbors = [];
    for (let dx = -1; dx <= 1; dx++) {
      for (let dy = -1; dy <= 1; dy++) {
        const key = this.hash(cx + dx, cy + dy);
        neighbors.push(...this.table[key]);
      }
    }
    return neighbors;
  }
}

Funkcja haszująca mapuje dowolne pary liczb całkowitych (cx, cy) na zakres [0, tableSize). Użycie dużych względnie pierwszych liczb zapobiega systematycznym kolizjom między sąsiednimi komórkami. Tablica musi być czyszczona i odbudowywana co klatkę, gdy cząstki się poruszają — czyszczenie całej tablicy zajmuje czas O(tableSize), więc tableSize powinno być wybrane jako niewielka wielokrotność n, a nie wykładniczo większe.

Projektowanie funkcji haszującej

Kolizje haszy (różne komórki mapujące się na to samo miejsce w tabeli) nie powodują błędów poprawności — sprawdzamy rzeczywiste odległości, by odfiltrować fałszywe trafienia — ale zwiększają nakład pracy. Słaba funkcja haszująca (np. (cx + cy) % tableSize) powoduje wiele kolizji między komórkami wzdłuż przekątnych. Iloczyn dużych liczb pierwszych rozprasza współrzędne komórek quasi-losowo w tabeli.

Dla shaderów obliczeniowych WebGL, gdzie arytmetyka na liczbach całkowitych jest ograniczona, powszechną alternatywą jest użycie pojedynczego dużego klucza całkowitego key = cx * MAX_Y + cy i przechowywanie cząstek w posortowanej tablicy według klucza, co umożliwia wyszukiwanie binarne. Pozwala to uniknąć dynamicznej alokacji pamięci (bez tablic JavaScript), co jest kluczowe dla pamięci GPU.

Wybór tableSize: przy n cząstkach i średniej liczbie sąsiadów k, każda komórka zawiera około k/9 cząstek. Ustawienie tableSize = 2n daje współczynnik wypełnienia 0,5 i akceptowalną liczbę kolizji. Mniejsze tabele ryzykują wiele kolizji (degradując do O(n²) w najgorszym przypadku); większe tabele marnują przepustowość pamięci podczas czyszczenia.

Podejście z sumą prefiksową dla równoległości GPU

Tablice haszujące JavaScript wykorzystują dynamiczne tablice i nie mogą działać równolegle na GPU. Implementacje WebGL wykorzystują podejście sumy prefiksowej (scan):

1. Przebieg zliczający: każda cząstka atomowo zwiększa licznik dla swojej komórki (shader obliczeniowy).
2. Suma prefiksowa: obliczenie sum skumulowanych, aby znaleźć indeks początkowy każdej komórki w posortowanej tablicy.
3. Przebieg wypełniający: zapis każdej cząstki do slotu jej komórki w posortowanej tablicy, wykorzystując sumę prefiksową jako wskaźnik zapisu.
4. Zapytanie: dla dowolnej komórki odczytaj cząstki od indeksu prefixSum[komórka] do prefixSum[komórka+1].

Tworzy to idealnie upakowaną, posortowaną tablicę cząstek według komórek — bez list wiązanych, bez dynamicznej alokacji, w pełni scalony dostęp do pamięci na GPU. Sam przebieg sumy prefiksowej działa w O(log n) równoległych krokach przy użyciu algorytmu równoległego skanu. Łącznie pełna odbudowa hasza przestrzennego kosztuje O(n) pracy, ale tylko O(log n) kroków równoległych — praktycznie za darmo w porównaniu do samych obliczeń sąsiedztwa.

Haszowanie przestrzenne kontra BVH: kiedy stosować które

• Haszowanie przestrzenne: najlepsze dla cząstek o jednakowym rozmiarze (jednolity rozmiar komórki jest optymalny), scen dynamicznych (odbudowa co klatkę to O(n)), dużej liczby cząstek (10 tys.–1 mln), rozkładów jednolitych. Stosowane w płynach SPH, piasku, stadach ptaków (boids).
• BVH: najlepsze dla obiektów o zmiennym rozmiarze (BVH dostosowuje się do rozmiarów obiektów), scen statycznych lub wolno zmieniających się (odbudowa jest kosztowna), zapytań o przecięcia promieni, scen z dużymi pustymi obszarami. Stosowane w ray tracingu, fizyce ciał sztywnych o złożonych kształtach.