Skip to content
GitLab
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
SWIFT
SWIFTsim
Commits
8d372771
Commit
8d372771
authored
Sep 26, 2016
by
Matthieu Schaller
Browse files
Written down the full definitions of the Minimal-SPH and Gadget2-SPH flavours.
parent
f104bb3a
Changes
5
Hide whitespace changes
Inline
Side-by-side
.gitignore
View file @
8d372771
...
...
@@ -83,6 +83,7 @@ theory/latex/swift.pdf
theory/SPH/kernel/kernels.pdf
theory/SPH/kernel/kernel_derivatives.pdf
theory/SPH/kernel/kernel_definitions.pdf
theory/SPH/flavours/sph_flavours.pdf
theory/paper_pasc/pasc_paper.pdf
m4/libtool.m4
...
...
theory/SPH/flavours/bibliography.bib
0 → 100644
View file @
8d372771
@ARTICLE
{
Price2012
,
author
=
{{Price}, D.~J.}
,
title
=
"{Smoothed particle hydrodynamics and magnetohydrodynamics}"
,
journal
=
{Journal of Computational Physics}
,
archivePrefix
=
"arXiv"
,
eprint
=
{1012.1885}
,
primaryClass
=
"astro-ph.IM"
,
year
=
2012
,
month
=
feb
,
volume
=
231
,
pages
=
{759-794}
,
doi
=
{10.1016/j.jcp.2010.12.011}
,
adsurl
=
{http://adsabs.harvard.edu/abs/2012JCoPh.231..759P}
,
adsnote
=
{Provided by the SAO/NASA Astrophysics Data System}
}
@ARTICLE
{
Springel2005
,
author
=
{{Springel}, V.}
,
title
=
"{The cosmological simulation code GADGET-2}"
,
journal
=
{\mnras}
,
eprint
=
{astro-ph/0505010}
,
keywords
=
{methods: numerical, galaxies: interactions, dark matter}
,
year
=
2005
,
month
=
dec
,
volume
=
364
,
pages
=
{1105-1134}
,
doi
=
{10.1111/j.1365-2966.2005.09655.x}
,
adsurl
=
{http://adsabs.harvard.edu/abs/2005MNRAS.364.1105S}
,
adsnote
=
{Provided by the SAO/NASA Astrophysics Data System}
}
@ARTICLE
{
Hopkins2013
,
author
=
{{Hopkins}, P.~F.}
,
title
=
"{A general class of Lagrangian smoothed particle hydrodynamics methods and implications for fluid mixing problems}"
,
journal
=
{\mnras}
,
archivePrefix
=
"arXiv"
,
eprint
=
{1206.5006}
,
primaryClass
=
"astro-ph.IM"
,
keywords
=
{hydrodynamics, instabilities, turbulence, methods: numerical, cosmology: theory}
,
year
=
2013
,
month
=
feb
,
volume
=
428
,
pages
=
{2840-2856}
,
doi
=
{10.1093/mnras/sts210}
,
adsurl
=
{http://adsabs.harvard.edu/abs/2013MNRAS.428.2840H}
,
adsnote
=
{Provided by the SAO/NASA Astrophysics Data System}
}
@ARTICLE
{
Springel2002
,
author
=
{{Springel}, V. and {Hernquist}, L.}
,
title
=
"{Cosmological smoothed particle hydrodynamics simulations: the entropy equation}"
,
journal
=
{\mnras}
,
eprint
=
{astro-ph/0111016}
,
keywords
=
{methods: numerical, galaxies: evolution, galaxies: starburst, methods: numerical, galaxies: evolution, galaxies: starburst}
,
year
=
2002
,
month
=
jul
,
volume
=
333
,
pages
=
{649-664}
,
doi
=
{10.1046/j.1365-8711.2002.05445.x}
,
adsurl
=
{http://adsabs.harvard.edu/abs/2002MNRAS.333..649S}
,
adsnote
=
{Provided by the SAO/NASA Astrophysics Data System}
}
@ARTICLE
{
Balsara1995
,
author
=
{{Balsara}, D.~S.}
,
title
=
"{von Neumann stability analysis of smooth particle hydrodynamics--suggestions for optimal algorithms}"
,
journal
=
{Journal of Computational Physics}
,
year
=
1995
,
volume
=
121
,
pages
=
{357-372}
,
doi
=
{10.1016/S0021-9991(95)90221-X}
,
adsurl
=
{http://adsabs.harvard.edu/abs/1995JCoPh.121..357B}
,
adsnote
=
{Provided by the SAO/NASA Astrophysics Data System}
}
theory/SPH/flavours/run.sh
0 → 100755
View file @
8d372771
#!/bin/bash
pdflatex sph_flavours.tex
bibtex sph_flavours.aux
pdflatex sph_flavours.tex
pdflatex sph_flavours.tex
theory/SPH/flavours/sph_flavours.tex
0 → 100644
View file @
8d372771
\documentclass
[fleqn, usenatbib, useAMS, a4paper]
{
mnras
}
\usepackage
{
graphicx
}
\usepackage
{
amsmath,paralist,xcolor,xspace,amssymb
}
\usepackage
{
times
}
\usepackage
[super]
{
nth
}
\newcommand
{
\swift
}{{
\sc
swift
}
\xspace
}
\newcommand
{
\gadget
}{{
\sc
gadget
}
\xspace
}
\newcommand
{
\dd
}
[2]
{
\frac
{
\partial
#1
}{
\partial
#2
}}
\renewcommand
{
\vec
}
[1]
{{
\mathbf
{
#1
}}}
\newcommand
{
\Wij
}{
\overline
{
\nabla
_
xW
_{
ij
}}}
%\setlength{\mathindent}{0pt}
%opening
\title
{
SPH flavours in SWIFT
}
\author
{
Matthieu Schaller
}
\begin{document}
\maketitle
For vector quantities, we define
$
\vec
{
a
}_{
ij
}
\equiv
(
\vec
{
a
}_
i
-
\vec
{
a
}_
j
)
$
. Barred quantities (
$
\bar
a
_{
ij
}$
) correspond to the
average over two particles of a quantity:
$
\bar
a
_{
ij
}
\equiv
\frac
{
1
}{
2
}
(
a
_
i
+
a
_
j
)
$
. To simplify notations, we also define the
vector quantity
$
\Wij
\equiv
\frac
{
1
}{
2
}
\left
(
W
(
\vec
{
x
}_{
ij
}
, h
_
i
)
+
\nabla
_
x W
(
\vec
{
x
}_{
ij
}
,h
_
j
)
\right
)
$
.
\subsection
{
Minimal SPH
}
\label
{
sec:sph:minimal
}
This is the simplest fully-conservative version of SPH using the
internal energy
$
u
$
as a thermal variable that can be written
down. Its implementation in
\swift
should be understood as a text-book
example and template for more advanced implementation. A full
derivation of the equations can be found in the review of
\cite
{
Price2012
}
. Our implementation follows their equations (27), (43),
(44), (45), (101), (103) and (104) with
$
\beta
=
3
$
and
$
\alpha
_
u
=
0
$
. We
summarize the equations here.
\subsubsection
{
Density and other fluid properties (
\nth
{
1
}
neighbour loop)
}
For a set of particles
$
i
$
with positions
$
\vec
{
x
}_
i
$
with velocities
$
\vec
{
v
}_
i
$
, masses
$
m
_
i
$
, sthermal energy per unit mass
$
u
_
i
$
and
smoothing length
$
h
_
i
$
, we compute the density for each particle:
\begin{equation}
\rho
_
i
\equiv
\rho
(
\vec
{
x
}_
i) =
\sum
_
j m
_
j W(
\vec
{
x
}_{
ij
}
, h
_
i),
\label
{
eq:sph:minimal:rho
}
\end{equation}
and the derivative of its density with respect to
$
h
$
:
\begin{equation}
\label
{
eq:sph:minimal:rho
_
dh
}
\rho
_{
\partial
h
_
i
}
\equiv
\dd
{
\rho
}{
h
}
(
\vec
{
x
}_
i) =
\sum
_
j m
_
j
\dd
{
W
}{
h
}
(
\vec
{
x
}_{
ij
}
, h
_
i).
\end{equation}
The gradient terms (``h-terms'') can then be computed from the density
and its derivative:
\begin{equation}
f
_
i
\equiv
\left
(1 +
\frac
{
h
_
i
}{
3
\rho
_
i
}
\rho
_{
\partial
h
_
i
}
\right
)
^{
-1
}
.
\label
{
eq:sph:minimal:f
_
i
}
\end{equation}
Using the pre-defined equation of state, the pressure
$
P
_
i
$
and the sound
speed
$
c
_
i
$
at the location of particle
$
i
$
can now be computed from
$
\rho
_
i
$
and
$
u
_
i
$
:
\begin{align}
P
_
i
&
= P
_{
\rm
eos
}
(
\rho
_
i, u
_
i),
\label
{
eq:sph:minimal:P
}
\\
c
_
i
&
= c
_{
\rm
eos
}
(
\rho
_
i, u
_
i).
\label
{
eq:sph:minimal:c
}
\end{align}
\subsubsection
{
Hydrodynamical accelerations (
\nth
{
2
}
neighbour loop)
}
We can then proceed with the second loop over
neighbours. The signal velocity
$
v
_{{
\rm
sig
}
,ij
}$
between two particles is given by
\begin{align}
\mu
_{
ij
}
&
=
\begin{cases}
\frac
{
\vec
{
v
}_{
ij
}
\cdot
\vec
{
x
}_{
ij
}}{
|
\vec
{
x
}_{
ij
}
|
}
&
\rm
{
if
}
~
\vec
{
v
}_{
ij
}
\cdot
\vec
{
x
}_{
ij
}
< 0,
\\
0
&
\rm
{
otherwise
}
,
\\
\end{cases}
\nonumber\\
v
_{{
\rm
sig
}
,ij
}
&
= c
_
i + c
_
j - 3
\mu
_{
ij
}
.
\label
{
eq:sph:minimal:v
_
sig
}
\end{align}
We also use these two quantities for the calculation of the viscosity term:
\begin{equation}
\nu
_{
ij
}
= -
\frac
{
1
}{
2
}
\frac
{
\alpha
\mu
_{
ij
}
v
_{{
\rm
sig
}
,ij
}}{
\bar\rho
_{
ij
}}
\label
{
eq:sph:minimal:nu
_
ij
}
\end{equation}
The fluid accelerations are then given by
\begin{align}
\frac
{
d
\vec
{
v
}_
i
}{
dt
}
= -
\sum
_
j m
_
j
&
\left
[
\frac
{
f
_
iP
_
i
}{
\rho
_
i
^
2
}
\nabla
_
x W(
\vec
{
x
}_{
ij
}
, h
_
i)
\nonumber
+
\frac
{
f
_
jP
_
j
}{
\rho
_
j
^
2
}
\nabla
_
x W(
\vec
{
x
}_{
ij
}
,h
_
j)
\right
.
\\
&
+
\left
.
\bigg
.
\nu
_{
ij
}
\Wij
\right
],
\label
{
eq:sph:minimal:dv
_
dt
}
\end{align}
and the change in internal energy,
\begin{align}
\frac
{
du
_
i
}{
dt
}
=
\sum
_
j m
_
j
&
\left
[
\frac
{
f
_
iP
_
i
}{
\rho
_
i
^
2
}
\vec
{
v
}_{
ij
}
\cdot
\nabla
_
x W(
\vec
{
x
}_{
ij
}
, h
_
i)
\right
.
\label
{
eq:sph:minimal:du
_
dt
}
\\
&
+
\left
.
\frac
{
1
}{
2
}
\nu
_{
ij
}
\vec
{
v
}_{
ij
}
\cdot\Big
.
\Wij\right
],
\nonumber
\end{align}
where in both cases the first line corresponds to the standard SPH
term and the second line to the viscuous accelerations.
We also compute an estimator of the change in smoothing length to be
used in the prediction step. This is an estimate of the local
divergence of the velocity field compatible with the accelerations
computed above:
\begin{equation}
\frac
{
dh
_
i
}{
dt
}
= -
\frac
{
1
}{
3
}
h
_
i
\sum
_
j
\frac
{
m
_
j
}{
\rho
_
j
}
\vec
{
v
}_{
ij
}
\cdot
\nabla
_
x W(
\vec
{
x
}_{
ij
}
, h
_
i).
\label
{
eq:sph:minimal:dh
_
dt
}
\end{equation}
and update the signal velocity of the particles:
\begin{equation}
v
_{{
\rm
sig
}
,i
}
=
\max
_
j
\left
( v
_{{
\rm
sig
}
,ij
}
\right
).
\label
{
eq:sph:minimal:v
_
sig
_
update
}
\end{equation}
All the quantities required for time integration have now been obtained.
\subsubsection
{
Time integration
}
For each particle, we compute a time-step given by the CFL condition:
\begin{equation}
\Delta
t = 2 C
_{
\rm
CFL
}
\frac
{
H
_
i
}{
v
_{{
\rm
sig
}
,i
}}
,
\label
{
eq:sph:minimal:dt
}
\end{equation}
where
$
C
_{
\rm
CFL
}$
is a free dimensionless parameter and
$
H
_
i
=
\gamma
h
_
i
$
is the
kernel support size. Particles can then be ``kicked'' forward in time:
\begin{align}
\vec
{
v
}_
i
&
\rightarrow
\vec
{
v
}_
i +
\frac
{
d
\vec
{
v
}_
i
}{
dt
}
\Delta
t
\label
{
eq:sph:minimal:kick
_
v
}
\\
u
_
i
&
\rightarrow
u
_
i +
\frac
{
du
_
i
}{
dt
}
\Delta
t
\label
{
eq:sph:minimal:kick
_
u
}
\\
P
_
i
&
\rightarrow
P
_{
\rm
eos
}
\left
(
\rho
_
i, u
_
i
\right
)
\label
{
eq:sph:minimal:kick
_
P
}
,
\\
c
_
i
&
\rightarrow
c
_{
\rm
eos
}
\left
(
\rho
_
i, u
_
i
\right
)
\label
{
eq:sph:minimal:kick
_
c
}
,
\end{align}
where we used the pre-defined equation of state to compute the new
value of the pressure and sound-speed.
\subsubsection
{
Particle properties prediction
}
Inactive particles need to have their quantities predicted forward in
time in the ``drift'' operation. We update them as follows:
\begin{align}
h
_
i
&
\rightarrow
h
_
i
\exp\left
(
\frac
{
1
}{
h
_
i
}
\frac
{
dh
_
i
}{
dt
}
\Delta
t
\right
),
\label
{
eq:sph:minimal:drift
_
h
}
\\
\rho
_
i
&
\rightarrow
\rho
_
i
\exp\left
(-
\frac
{
3
}{
h
_
i
}
\frac
{
dh
_
i
}{
dt
}
\Delta
t
\right
),
\label
{
eq:sph:minimal:drift
_
rho
}
\\
P
_
i
&
\rightarrow
P
_{
\rm
eos
}
\left
(
\rho
_
i, u
_
i +
\frac
{
du
_
i
}{
dt
}
\Delta
t
\right
),
\label
{
eq:sph:minimal:drift
_
P
}
\\
c
_
i
&
\rightarrow
c
_{
\rm
eos
}
\left
(
\rho
_
i, u
_
i +
\frac
{
du
_
i
}{
dt
}
\Delta
t
\right
)
\label
{
eq:sph:minimal:drift
_
c
}
,
\end{align}
where, as above, the last two updated quantities are obtained using
the pre-defined equation of state. Note that the thermal energy
$
u
_
i
$
itself is
\emph
{
not
}
updated.
\subsection
{
Gadget-2 SPH
}
\label
{
sec:sph:gadget2
}
This flavour of SPH is the one implemented in the
\gadget
-2 code
\citep
{
Springel2005
}
. The basic equations were derived by
\cite
{
Springel2002
}
and also includes a
\cite
{
Balsara1995
}
switch for
the suppression of viscosity. The implementation here follows closely the
presentation of
\cite
{
Springel2005
}
. Specifically, we use their equations (5), (7),
(8), (9), (10), (13), (14) and (17). We summarize them here for completeness.
\subsubsection
{
Density and other fluid properties (
\nth
{
1
}
neighbour loop)
}
For a set of particles
$
i
$
with positions
$
\vec
{
x
}_
i
$
with velocities
$
\vec
{
v
}_
i
$
, masses
$
m
_
i
$
, entropic function per unit mass
$
A
_
i
$
and
smoothing length
$
h
_
i
$
, we compute the density, derivative of the density with respect
to
$
h
$
and the ``f-terms'' in a similar way to the minimal-SPH case
(Eq.
\ref
{
eq:sph:minimal:rho
}
,
\ref
{
eq:sph:minimal:rho
_
dh
}
and
\ref
{
eq:sph:minimal:f
_
i
}
). From these the pressure and sound-speed can
be computed using the predefined equation of state:
\begin{align}
P
_
i
&
= P
_{
\rm
eos
}
(
\rho
_
i, A
_
i),
\label
{
eq:sph:gadget2:P
}
\\
c
_
i
&
= c
_{
\rm
eos
}
(
\rho
_
i, A
_
i).
\label
{
eq:sph:gadget2:c
}
\end{align}
We additionally compute the divergence and
curl of the velocity field using standard SPH expressions:
\begin{align}
\nabla\cdot\vec
{
v
}_
i
\equiv\nabla\cdot
\vec
{
v
}
(
\vec
{
x
}_
i)
&
=
\frac
{
1
}{
\rho
_
i
}
\sum
_
j m
_
j
\vec
{
v
}_{
ij
}
\cdot\nabla
_
x W(
\vec
{
x
}_{
ij
}
, h
_
i)
\label
{
eq:sph:gadget2:div
_
v
}
,
\\
\nabla\times\vec
{
v
}_
i
\equiv
\nabla\times
\vec
{
v
}
(
\vec
{
x
}_
i)
&
=
\frac
{
1
}{
\rho
_
i
}
\sum
_
j m
_
j
\vec
{
v
}_{
ij
}
\times\nabla
_
x W(
\vec
{
x
}_{
ij
}
, h
_
i)
\label
{
eq:sph:gadget2:rot
_
v
}
.
\end{align}
These are used to construct the
\cite
{
Balsara1995
}
switch
$
B
_
i
$
:
\begin{equation}
B
_
i =
\frac
{
|
\nabla\cdot\vec
{
v
}_
i|
}{
|
\nabla\cdot\vec
{
v
}_
i| +
|
\nabla\times\vec
{
v
}_
i| + 10
^{
-4
}
c
_
i / h
_
i
}
,
\label
{
eq:sph:gadget2:balsara
}
\end{equation}
where the last term in the denominator is added to prevent numerical instabilities.
\subsubsection
{
Hydrodynamical accelerations (
\nth
{
2
}
neighbour loop)
}
The accelerations are computed in a similar fashion to the minimal-SPH
case with the exception of the viscosity term which is now modified to
include the switch. Instead of Eq.
\ref
{
eq:sph:minimal:nu
_
ij
}
, we get:
\begin{equation}
\nu
_{
ij
}
= -
\frac
{
1
}{
2
}
\frac
{
\alpha
\bar
B
_{
ij
}
\mu
_{
ij
}
v
_{{
\rm
sig
}
,ij
}}{
\bar\rho
_{
ij
}}
,
\label
{
eq:sph:gadget2:nu
_
ij
}
\end{equation}
whilst equations
\ref
{
eq:sph:minimal:v
_
sig
}
,
\ref
{
eq:sph:minimal:dv
_
dt
}
,
\ref
{
eq:sph:minimal:dh
_
dt
}
and
\ref
{
eq:sph:minimal:v
_
sig
_
update
}
remain unchanged. The only other
change is the equation of motion for the thermodynamic variable which
now has to be describing the evolution of the entropic function and
not the evolution of the thermal energy. It reads:
\begin{equation}
\frac
{
dA
_
i
}{
dt
}
=
\frac
{
1
}{
2
}
\frac
{
\gamma
-1
}{
\rho
^{
\gamma
-1
}}
\sum
_
j
m
_
j
\nu
_{
ij
}
\vec
{
v
}_{
ij
}
\cdot
\Wij
\end{equation}
\subsubsection
{
Time integration
}
The time-step condition is identical to the Minimal-SPH case
(Eq.
\ref
{
eq:sph:minimal:dt
}
). The same applies to the integration
forward in time (Eq.
\ref
{
eq:sph:minimal:kick
_
v
}
to
\ref
{
eq:sph:minimal:kick
_
c
}
) with the exception of the change in
internal energy (Eq.
\ref
{
eq:sph:minimal:kick
_
u
}
) which gets replaced
by an integration for the the entropy:
\begin{align}
\vec
{
v
}_
i
&
\rightarrow
\vec
{
v
}_
i +
\frac
{
d
\vec
{
v
}_
i
}{
dt
}
\Delta
t
\label
{
eq:sph:gadget2:kick
_
v
}
\\
A
_
i
&
\rightarrow
A
_
i +
\frac
{
dA
_
i
}{
dt
}
\Delta
t
\label
{
eq:sph:gadget2:kick
_
A
}
\\
P
_
i
&
\rightarrow
P
_{
\rm
eos
}
\left
(
\rho
_
i, A
_
i
\right
)
\label
{
eq:sph:gadget2:kick
_
P
}
,
\\
c
_
i
&
\rightarrow
c
_{
\rm
eos
}
\left
(
\rho
_
i, A
_
i
\right
)
\label
{
eq:sph:gadget2:kick
_
c
}
,
\end{align}
where, once again, we made use of the equation of state relating
thermodynaimcal quantities.
\subsubsection
{
Particle properties prediction
}
The prediction step is also identical to the Minimal-SPH case with the
entropic function replacing the thermal energy.
\begin{align}
h
_
i
&
\rightarrow
h
_
i
\exp\left
(
\frac
{
1
}{
h
_
i
}
\frac
{
dh
_
i
}{
dt
}
\Delta
t
\right
),
\label
{
eq:sph:gadget2:drift
_
h
}
\\
\rho
_
i
&
\rightarrow
\rho
_
i
\exp\left
(-
\frac
{
3
}{
h
_
i
}
\frac
{
dh
_
i
}{
dt
}
\Delta
t
\right
),
\label
{
eq:sph:gadget2:drift
_
rho
}
\\
P
_
i
&
\rightarrow
P
_{
\rm
eos
}
\left
(
\rho
_
i, A
_
i +
\frac
{
dA
_
i
}{
dt
}
\Delta
t
\right
),
\label
{
eq:sph:gadget2:drift
_
P
}
\\
c
_
i
&
\rightarrow
c
_{
\rm
eos
}
\left
(
\rho
_
i, A
_
i +
\frac
{
dA
_
i
}{
dt
}
\Delta
t
\right
)
\label
{
eq:sph:gadget2:drift
_
c
}
,
\end{align}
where, as above, the last two updated quantities are obtained using
the pre-defined equation of state. Note that the entropic function
$
A
_
i
$
itself is
\emph
{
not
}
updated.
\subsection
{
Pressure-Entropy SPH
}
\cite
{
Hopkins2013
}
\label
{
sec:sph:pe
}
\subsubsection
{
Hydrodynamical accelerations
}
\subsubsection
{
Particle properties prediction
}
\subsection
{
Pressure-Energy SPH
}
\label
{
sec:sph:pu
}
\cite
{
Hopkins2013
}
\\
to be done
\subsection
{
Anarchy SPH
}
\label
{
sec:sph:anarchy
}
to be done
\bibliographystyle
{
mnras
}
\bibliography
{
./bibliography.bib
}
\end{document}
theory/SPH/kernel/kernel_definitions.tex
View file @
8d372771
\documentclass
[usenatbib, useAMS,a4paper]
{
mnras
}
\documentclass
[
fleqn,
usenatbib, useAMS,a4paper]
{
mnras
}
\usepackage
{
graphicx
}
\usepackage
{
amsmath,paralist,xcolor,xspace,amssymb
}
\usepackage
{
times
}
...
...
@@ -27,14 +27,12 @@ The desirable properties of an SPH kernels $W(\vec{x},h)$ are:
\item
$
W
(
\vec
{
x
}
,h
)
$
should have a finite support and be cheap to
compute.
\end{enumerate}
As a consequence, the smoothing kernels used in SPH can
hence be written (in 3D) as
\begin{equation}
W(
\vec
{
x
}
,h)
\equiv
\frac
{
1
}{
H
^
3
}
f
\left
(
\frac
{
|
\vec
{
x
}
|
}{
H
}
\right
),
\end{equation}
where
$
H
=
\gamma
h
$
is defined below and
$
f
(
u
)
$
is a dimensionless
function, usually a low-order polynomial, such that
$
f
(
u
\geq
1
)
=
0
$
and normalised such that
...
...
@@ -42,7 +40,6 @@ and normalised such that
\begin{equation}
\int
f(|
\vec
{
u
}
|)
{
\rm
d
}^
3u = 1.
\end{equation}
$
H
$
is the kernel's support radius and is used as the ``smoothing
length'' in the Gadget code(
{
i.e.
}
$
H
=
h
$
). This definition is,
however, not very physical and makes comparison of kernels at a
...
...
@@ -55,13 +52,11 @@ standard deviation
\sigma
^
2
\equiv
\frac
{
1
}{
3
}
\int
\vec
{
u
}^
2 W(
\vec
{
u
}
,h)
{
\rm
d
}^
3u.
\label
{
eq:sph:sigma
}
\end{equation}
The smoothing length is then:
\begin{equation}
h
\equiv
2
\sigma
.
\label
{
eq:sph:h
}
\end{equation}
Each kernel,
{
\it
i.e.
}
defintion of
$
f
(
u
)
$
, will have a different
ratio
$
\gamma
=
H
/
h
$
. So for a
\emph
{
fixed resolution
}
$
h
$
, one will
have different kernel support sizes,
$
H
$
, and a different number of
...
...
@@ -72,18 +67,15 @@ separation:
\begin{equation}
h =
\eta
\langle
x
\rangle
=
\eta
\left
(
\frac
{
m
}{
\rho
}
\right
)
^{
1/3
}
,
\end{equation}
where
$
\rho
$
is the local density of the fluid and
$
m
$
the SPH
particle mass.
The (weighted) number of neighbours within the kernel support is a
useful quantity to use in implementations of SPH. It is defined
as
(in
3D)
useful quantity to use in implementations of SPH. It is defined (in
3D)
as:
\begin{equation}
N
_{
\rm
ngb
}
\equiv
\frac
{
4
}{
3
}
\pi
\left
(
\frac
{
H
}{
h
}
\eta\right
)
^
3.
\end{equation}
Once the fixed ratio
$
\gamma
=
H
/
h
$
is known (via equations
\ref
{
eq:sph:sigma
}
and
\ref
{
eq:sph:h
}
) for a given kernel, the number
of neighbours only depends on the resolution parameter
$
\eta
$
. For
...
...
@@ -96,23 +88,22 @@ resolution to $\eta=1.2348$ yields the commonly used value $N_{\rm
The
\swift
kernels are split into two categories, the B-splines
(
$
M
_{
4
,
5
,
6
}$
) and the Wendland kernels (
$
C
2
$
,
$
C
4
$
and
$
C
6
$
). In all
cases we impose
$
f
(
u>
1
)
=
0
$
.
\\
The spline kernels are defined as:
\begin{align}
f(u)
&
= C M
_
n(u),
\\
f(u)
&
= C M
_
n(u),
\nonumber
\\
M
_
n(u)
&
\equiv
\frac
{
1
}{
2
\pi
}
\int
_{
-
\infty
}^{
\infty
}
\left
(
\frac
{
\sin\left
(k/n
\right
)
}{
k/n
}
\right
)
^
n
\cos\left
(ku
\right
)
{
\rm
d
}
k,
d
}
k
\nonumber
,
\end{align}
whilst the Wendland kernels read
whilst the Wendland kernels are constructed from:
\begin{align}
f(u)
&
= C
\Psi
_{
i,j
}
(u),
\\
\Psi
_{
i,j
}
(u)
&
\equiv
\mathcal
{
I
}^
k
\left
[\max\left(1-u,0\right)\right]
,
\\
\mathcal
{
I
}
[f](u)
&
\equiv
\int
_
u
^
\infty
f
\left
(k
\right
)k
{
\rm
d
}
k.
f(u)
&
= C
\Psi
_{
i,j
}
(u),
\nonumber\\
\Psi
_{
i,j
}
(u)
&
\equiv
\mathcal
{
I
}^
k
\left
[\max\left(1-u,0\right)\right]
,
\nonumber\\
\mathcal
{
I
}
[f](u)
&
\equiv
\int
_
u
^
\infty
f
\left
(k
\right
)k
{
\rm
d
}
k.
\nonumber
\end{align}
\subsubsection
{
Cubic spline (
$
M
_
4
$
) kernel
}
...
...
@@ -230,11 +221,11 @@ All kernels available in \swift are shown on Fig.~\ref{fig:sph:kernels}.
The derivatives of the kernel function have relatively simple
expressions and are shown on Fig.~
\ref
{
fig:sph:kernel
_
derivatives
}
:
\begin{
eqnarray*
}
\vec\nabla
_
x W(
\vec
{
x
}
,h)
&
=
&
\frac
{
1
}{
h
^
4
}
f'
\left
(
\frac
{
|
\vec
{
x
}
|
}{
h
}
\right
)
\frac
{
\vec
{
x
}}{
|
\vec
{
x
}
|
}
,
\\
\frac
{
\partial
W(
\vec
{
x
}
,h)
}{
\partial
h
}
&
=
&
-
\frac
{
1
}{
h
^
4
}
\left
[3f
\left
(
\frac
{
|
\vec
{
x
}
|
}{
h
}
\right
) +
\begin{
align
}
\vec\nabla
_
x W(
\vec
{
x
}
,h)
&
=
\frac
{
1
}{
h
^
4
}
f'
\left
(
\frac
{
|
\vec
{
x
}
|
}{
h
}
\right
)
\frac
{
\vec
{
x
}}{
|
\vec
{
x
}
|
}
,
\\
\frac
{
\partial
W(
\vec
{
x
}
,h)
}{
\partial
h
}
&
=-
\frac
{
1
}{
h
^
4
}
\left
[3f
\left
(
\frac
{
|
\vec
{
x
}
|
}{
h
}
\right
) +
\frac
{
|
\vec
{
x
}
|
}{
h
}
f'
\left
(
\frac
{
|
\vec
{
x
}
|
}{
h
}
\right
)
\right
].
\end{
eqnarray*
}
\end{
align
}
Note that for all the kernels listed above,
$
f'
(
0
)
=
0
$
.
\bibliographystyle
{
mnras
}
...
...
Write
Preview
Supports
Markdown
0%
Try again
or
attach a new file
.
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment