カーネル関数 - YuriOku/1D_SPH GitHub Wiki

SPH法ではカーネル近似によって場を表現します。カーネル近似の下では、物理量 $f$ は

$f(x) = \int f(x') W(|x - x'|, h) dx'$

と表されます。ここで $W$ は粒子の広がりを表すカーネル関数、 $h$ はカーネル関数の広がりを表すスムージング長です。

カーネル関数は以下の3つの性質を満たさなければいけません：

規格化されている $\int W(\vec{r}, h) dV = 1$
$h \rightarrow 0$ の極限で $\delta$ 関数になる。
2回微分可能

カーネル近似

SPH粒子の半径にあたる長さ(カーネル関数の値が0になる距離)はカーネルサポート半径と呼ばれ、一般的にはカーネルサポート半径 $H = 2h$ です。ただし、宇宙物理学の分野でよく使われているSPHシミュレーションコードGADGET (Springel et al. 2001; Springel 2005; Springel et al. 2020)では $H = h$ が採用されており、1D_SPHでもこの定義を採用しています。通常3次元シミュレーションでは、カーネルサポート半径内に含まれる粒子の数が一定値 $N_H$ となるようにスムージング長を決めます。 1次元コードである1D_SPHでは、

$h = \eta \left(\frac{m}{\rho}\right)$

となる $h$ を反復法によって求めます。ここで $\eta$ はパラメータで、標準値は $\eta = 2.4$ としています。

Cubic spline 関数

最も標準的なカーネル関数は3次スプライン関数 (Monaghan & Lattanzio 1985)

$W(r, h) = C \begin{cases} 1 - 6\left(\frac{r}{h}\right)^2 + 6\left(\frac{r}{h}\right)^3 & \left( 0 < \frac{r}{h} \leq \frac{1}{2} \right)\ 2\left( 1- \frac{r}{h} \right)^3 &\left( \frac{1}{2} < \frac{r}{h} \leq 1\right)\ 0 &\left( 1 < \frac{r}{h} \right) \end{cases}$

です。ここで $C$ は規格化定数で、1次元では $C = \frac{4}{3h}$ です。カーネル関数の半径微分は、

$\frac{\partial W(r, h)}{\partial r} = \frac{C}{h} \begin{cases} -12\left(\frac{r}{h}\right) + 18\left(\frac{r}{h}\right)^2 & \left( 0 < \frac{r}{h} \leq \frac{1}{2} \right)\ 6\left( 1- \frac{r}{h} \right)^2 &\left( \frac{1}{2} < \frac{r}{h} \leq 1\right)\ 0 &\left( 1 < \frac{r}{h} \right) \end{cases}$

です。スプライン関数シリーズには、より高次な Quartic spline (4次), Quintic spline (5次),... 関数もあります。

Wendland C2 関数

Wendland C2関数 (Wendland 1995)

$W(r, h) = C \begin{cases} \left(1 - \frac{r}{h} \right)^3 \left( 1 + 3\frac{r}{h} \right) & \left( \frac{r}{h} \leq 1 \right) \ 0 & \left( \frac{r}{h} > 1 \right) \end{cases}$

もよく使われるカーネル関数です。1次元では規格化定数は $C = \frac{5}{4h}$ です。 Wendland カーネルは、spline カーネルで近傍粒子数を大きくしたときに見られる pairing instability と呼ばれる数値不安定性を避けることができます (Dehnen & Aly 2012)。 Wendland C2関数の微分は

$\frac{\partial W(r, h)}{\partial r} = \frac{C}{h} \begin{cases} -\frac{4r}{h}\left(1 - \frac{r}{h} \right)^2 & \left( \frac{r}{h} \leq 1 \right) \ 0 & \left( \frac{r}{h} > 1 \right) \end{cases}$

です。Wendland シリーズにも、より高次なC4, C6関数があります。

参考文献

Dehnen, W. and Aly, H., “Improving convergence in smoothed particle hydrodynamics simulations without pairing instability”, Monthly Notices of the Royal Astronomical Society, vol. 425, no. 2, pp. 1068–1082, 2012. doi:10.1111/j.1365-2966.2012.21439.x.
Monaghan, J. J. and Lattanzio, J. C., “A refined particle method for astrophysical problems”, Astronomy and Astrophysics, vol. 149, no. 1, pp. 135–143, 1985.
Springel, V., Yoshida, N., and White, S. D. M., “GADGET: a code for collisionless and gasdynamical cosmological simulations”, New Astronomy, vol. 6, no. 2, pp. 79–117, 2001. doi:10.1016/S1384-1076(01)00042-2.
Springel, V., “The cosmological simulation code GADGET-2”, Monthly Notices of the Royal Astronomical Society, vol. 364, no. 4, pp. 1105–1134, 2005. doi:10.1111/j.1365-2966.2005.09655.x.
Springel, V., Pakmor, R., Zier, O., and Reinecke, M., “Simulating cosmic structure formation with the GADGET-4 code”, arXiv e-prints, 2020.
Wendland, H. "Piecewise polynomial, positive definite and compactly supported radial functions of minimal degree." Advances in computational Mathematics 4.1, 389-396, 1995