Metody numeryczne — jak fizyka staje się kodem

Problem: komputery nie potrafią liczyć rachunku różniczkowego

Równania fizyki opisują zmianę za pomocą pochodnych: prędkość to szybkość zmiany położenia, przyspieszenie to szybkość zmiany prędkości. W matematyce ciągłej to rachunek różniczkowy. W komputerze czas trzeba pokroić na małe dyskretne kroki (nazywane krokiem czasowym, Δt). Pytanie brzmi: jak dokładnie możemy przybliżyć ciągły ruch tymi drobnymi krokami?

Odpowiedź zależy od metody całkowania. Słaba metoda dryfuje w czasie — energia sztucznie narasta, aż symulowane sprężyny eksplodują. Dobra metoda utrzymuje całkowitą energię układu niemal stałą przez tysiące kroków.

Kluczowa obserwacja: każda symulacja fizyczna na tej stronie uruchamia jedną z tych trzech metod (lub ich kombinację) setki razy na sekundę. Wybór bezpośrednio wpływa na to, czy wahadło pozostanie stabilne przez 10 minut, czy zacznie oscylować w nieskończoność po 30 sekundach.

Trzy metody

Całkowanie Eulera

⚡ Najszybsza

Najprostsza metoda: użyj bieżącej prędkości, aby oszacować, gdzie znajdziesz się po Δt, a następnie zaktualizuj stan. Jedna ewaluacja funkcji na krok. Szybka, ale kumuluje błąd — energia powoli rośnie, przez co orbity rozkręcają się na zewnątrz, a sprężyny eksplodują.

// Euler jawny
position += velocity * dt;
velocity += acceleration * dt;

Zastosowanie: iskry cząsteczkowe, efekty wizualne, gdzie dokładność nie ma znaczenia

Runge-Kutta 4 (RK4)

🎯 Najdokładniejsza

Wykonuje cztery oszacowania pochodnej w obrębie każdego kroku czasowego i liczy z nich średnią ważoną. Znacznie lepsza dokładność niż Euler kosztem zaledwie 4-krotnie większych obliczeń. Standardowy wybór dla balistyki, mechaniki orbitalnej i symulacji chemicznych.

// Uproszczone RK4 dla pozycji/prędkości
const k1 = f(t,  state);
const k2 = f(t + dt/2, state + k1*dt/2);
const k3 = f(t + dt/2, state + k2*dt/2);
const k4 = f(t + dt,   state + k3*dt);
state += (k1 + 2*k2 + 2*k3 + k4) * dt/6;

Zastosowanie: symulacje orbitalne, balistyka, dyfuzja leków, kolumny płynu

Całkowanie Verleta / Leapfrog

🔒 Najstabilniejsza

Sprytna sztuczka: prędkość i położenie są przesunięte o pół kroku czasowego. Dzięki temu całkowita energia układu pozostaje niemal stała nawet w bardzo długich symulacjach — idealna dla dynamiki molekularnej i tkanin, gdzie symulacja musi działać bez końca. Używana wewnętrznie przez Cannon-es (silnik fizyczny na tej stronie).

// Velocity Verlet
posNext = pos + vel*dt + 0.5*acc*dt*dt;
velNext = vel + 0.5*(acc + accNext)*dt;

Zastosowanie: tkaniny, liny, dynamika molekularna, fizyka ciał sztywnych

Które symulacje używają której metody?

Balistyka — RK4 z oporem atmosferycznym; Euler nagromadziłby wystarczająco dużo błędu, by chybić cel o metry.
Symulacja tkaniny — Velocity Verlet przez Cannon-es; tkanina eksplodowałaby pod długotrwałym naprężeniem przy użyciu Eulera.
Gwiazdy podwójne — Euler symplektyczny (wariant lepiej zachowujący energię niż Euler jawny) dla długoterminowej stabilności orbitalnej.
Dyfuzja leku — RK4 na układzie równań różniczkowych; przedziały stężeń muszą być śledzone dokładnie przez wiele godzin.
Boids — Euler jawny; siły sterujące i tak zmieniają się co klatkę, więc dokładność ponad Δt nie ma znaczenia.

Wzorzec stałego kroku czasowego

Jedna dodatkowa komplikacja: monitory odświeżają się z częstotliwością 60 Hz, ale klatki nie zawsze renderują się dokładnie w 1/60 sekundy. Jeśli krok fizyki równa się krokowi renderowania, wolna klatka oznacza duże Δt, co może spowodować eksplozję sztywnych sprężyn.

Rozwiązaniem — stosowanym we wszystkich złożonych symulacjach na tej stronie — jest stały krok czasowy z akumulacją:

const FIXED_DT = 1 / 120; // fizyka 120 Hz
let accumulator = 0;

function tick(realDt) {
  accumulator += realDt;
  while (accumulator >= FIXED_DT) {
    physics.step(FIXED_DT);
    accumulator -= FIXED_DT;
  }
  render(); // działa z częstotliwością wyświetlacza
}

Fizyka zawsze postępuje identycznymi, bezpiecznymi krokami, niezależnie od tego, jak długo trwa renderowanie. Render po prostu rysuje aktualny stan fizyki.