SPH

Z Wikipedii, wolnej encyklopedii
Skocz do: nawigacja, szukaj

SPH (Smoothed Particle Hydrodynamics - wygładzona hydrodynamika cząstek) – jest to metoda numeryczna służąca do przeprowadzania symulacji numerycznych zachowania się płynów.

Metoda ta została zaproponowana w roku 1977 niezależnie przez R.A. Gingolda & J.J. Monaghana oraz przez L.B.Lucy do prowadzenia obliczeń z dziedziny astrofizyki. Początkowo używana była do symulacji ściśliwych cieczy nielepkich. Z czasem została rozwinięta na nieściśliwe ciecze lepkie znajdujące się w polu grawitacyjnym, a także na zagadnienia z dziedziny magnetohydrodynamiki.


Równania ruchu[edytuj | edytuj kod]

W metodzie SPH do opisu stanu cieczy używa się opisu Lagrange'a, gdzie siatka obliczeniowa porusza się wraz z przepływem cieczy. W takim wypadku równanie Naviera-Stokesa dla i-tej cząstki przybiera postać

 \frac{d\vec v}{dt} =-\frac{\nabla p}{\rho}+\vec a^{visc}+\vec \Phi

Gdzie

\vec v -prędkość
p - ciśnienie
\rho - gęstość
\vec a^{visc} - przyspieszenie wynikające z istnienia sił lepkości
\vec \Phi - przyspieszenie wynikające z obecności sił masowych (np. pola grawitacyjnego)

Metoda SPH[edytuj | edytuj kod]

Podstawy[edytuj | edytuj kod]

Metoda ta opiera się na teorii interpolacji. Ciągłe rozkłady takich parametrów jak gęstość czy ciśnienie cieczy zastępuje się odpowiednimi estymatami przy założonym pewnym jądrze interpolacji. Obliczenia wykonujemy dla dyskretnego zbioru N cząstek płynu.

Estymata pewnej wielkości A w pozycji i-tej cząstki jest dana jako

A_i=\sum_{j=1}^{N}m_j\frac{A_j}{\rho_j}W_{ij}

Natomiast estymata gradientu wielkości A jako

\nabla_iA_i=\sum_{j=1}^{N}m_j\frac{A_j}{\rho_j}\nabla_iW_{ij}

Gdzie

\rho_j - gęstość j-tej cząstki
m_j - masa j-tej cząstki
W_{ij}=W(r,h) - jądro interpolacji

Ilość sąsiadów i długość wygładzania[edytuj | edytuj kod]

Parametr h jest nazywany długością wygładzania (smoothing length). Jest to wielkość, która określa na jaką odległość cząstka może oddziaływać z innymi cząstkami. Najczęściej w symulacji pozostaje stała podczas trwania obliczeń. Należy tylko uwzględnić, aby w promieniu 2h znajdowała się odpowiednia liczba sąsiadów. Liczba ta powinna się wahać w granicach od N_N/2 do 2N_N.

Tabela 1. Ilość sąsiadów w promieniu 2h w zależności od liczby wymiarów symulacji
Jeden wymiar Dwa wymiary Trzy wymiary
liczba sąsiadów N_N 5 15 55

Jądro interpolacji[edytuj | edytuj kod]

Funkcja jądra interpolacji powinna być wybrana w postaci

 W(r,h)=\frac{1}{h^\nu}f(u)

gdzie

\nu - liczba wymiarów
u=r/h
r - odległość między cząstkami

Dodatkowo funkcja f(u) powinna spełniać warunki

\int f(u)dV=1

\lim_{h \rightarrow 0}f(u)=\delta(r)

przy czym dV jest elementem objętości, równym odpowiednio du, 2\pi u du lub 4\pi u^2 du w jednym, dwóch lub trzech wymiarach.

Najczęściej stosuje sie jednak jądro interpolacyjne zaproponowane przez Monaghana

W(r,h)=\frac{\sigma}{h^\nu}\left \{
\begin{matrix} 
1-3/2u^2+3/4u^3 & \quad 0\le u\le 1\\ 
1/4(2-u)^3 & \quad 1<u \le 2\\
0 & \quad u > 2 
\end{matrix}
\right.

A jego gradient w postaci

W(r,h)=\frac{\sigma}{h^{\nu+1}}\hat r\left \{
\begin{matrix} 
-3u+9/4u^2 & \quad 0 \le u \le 1 \\
-3/4(2-u)^2 &\quad 1<u \le 2\\
0 &\quad u>2
\end{matrix}
\right.

gdzie

\hat r=\frac{\vec r_i-\vec r_j}{|\vec r_i-\vec r_j|}
Tabela 2. Wartości parametrów \sigma i \nu w zależności od liczby wymiarów symulacji
Jeden wymiar Dwa wymiary Trzy wymiary
\nu 1 2 3
\sigma \frac{2}{3} \frac{10}{7\pi} \frac{1}{\pi}

Sztuczna lepkość[edytuj | edytuj kod]

Lepkość w metodzie SPH jest uwzględniana poprzez dodanie w równaniu Naviera-Stokesa przyspieszenia w postaci

\vec a_i^{visc}=-\sum_{j=1}^N m_j \Pi_{ij}\nabla_iW_{ij}

gdzie

\Pi_{ij}=\frac{-\alpha c \mu_{ij}+\beta \mu_{ij}^2}{\overline {\varrho}_{ij}}
\alpha, \beta - stałe, często przyjmowane jako \alpha = 0.5 \beta = 1.0
\overline {\varrho}_{ij}=(\varrho_{i}+\varrho_{j})/2


\mu_{ij}=\left\{ 
\begin{matrix}
\frac{h^2\vec v_{ij} \cdot \vec r_{ij}}{r_{ij}^2+\eta^2h^2} & \vec v_{ij} \cdot \vec r_{ij} <0 \\
0 &\vec v_{ij} \cdot \vec r_{ij} \ge 0
\end{matrix}
\right.
\eta -stała, często przyjmowane jako \eta= 0.1
\vec v_{ij}=\vec v_i-\vec v_j
\vec r_{ij}=\vec r_i-\vec r_j
r_{ij}=|\vec r_i-\vec r_j|

Równanie stanu[edytuj | edytuj kod]

Najczęściej stosuje się równanie stanu w postaci

c=\sqrt{\frac{p \gamma}{\rho}}

gdzie

\gamma=c_p/c_V
c - prędkość dźwięku w tym opisywanym płynie

przy założeniu warunku

 v_i \ll c

Innym równaniem stanu stosowanym dla cieczy jest

p=p_0\Bigg[\Bigg( \frac{\rho}{\rho _0}\Bigg)^{\gamma} -1\Bigg]

gdzie

\rho_0 i p_0 - maksymalne przyjęte wartości gęstości i ciśnienia
\gamma=7

Ostateczne równanie[edytuj | edytuj kod]

Po podstawieniu estymat do równania Naviera - Stokesa (i nie uwzględniając pola grawitacyjnego) otrzymujemy do rozwiązania równanie dla i-tej cząstki

\vec a_i=\frac{d\vec v_i}{dt}=-\frac{c^2}{\gamma}\sum_{j=1}^N m_j\Bigg(\frac{1}{\rho_i}+\frac{1}{\rho_j}+\Pi_{ij}\Bigg)\nabla_iW_{ij}

które można wyznaczać dowolną metodą numerycznego całkowania równań różniczkowych, ale przy założeniu, że krok czasowy spełnia warunek Couranta

\Delta t \le 0.25\frac{h}{c}

Literatura[edytuj | edytuj kod]

J.J. Monaghan, "Notes Smoothed Particle Hydrodynamics"

B.Schlatter , "A pedagogical tool using Smoothed Particle Hydrodynamics to model fluid flow past a system of cylinders"

J.P. Morris, "Analysis of Smoothed Particle Hydrodynamics with Applications"

M. Ellero, "Smoothed Particle Dynamics Methods for the Simulation of Viscoelastic Fluids"