2.2 Hydrodynamics
3.1.8 Equations of state
av
=−X
j
mjΠij∇iWij, (3.52) whereWij = (Wij(hi) +Wij(hj))/2. The corresponding heat production rate that arises from AV is also similar to that of the physical viscosity in Eq. (3.50) and given by
dui
dt av
= 1 2
X
j
mjΠijvij · ∇iWij. (3.53) Although this method is strongly recommended to apply to compressional fluid motion, spurious unexpected forces occur in shear flows [29]. Therefore, it is necessary to significantly reduce the AV in undesired regions, where one possible way is multiplication of Πij with the factorfij = (fi+fj)/2, which is the so called
“Balsara-switch” [56]. The factorsfi are given by fi = |∇ ·v|i
|∇ ·v|i+|∇ ×v|i+εi
, (3.54)
where εi = 10−4ci/hi 1 is due to numerical reasons, preventing a singularity.
In this formulation the switch reasonably suppresses artificial viscosity in regions of pure shear whereas the factor fi≈1 in the absence of shear flow.
The AV term produces heat from kinetic energy. However, it is also possible that discontinuities occur only across internal energy, which can not be smoothed by the AV procedure. For this purpose, an artificial conduction (AC) term can be constructed in analogy to the physical thermal conduction Eq. (3.46) (see e.g.
Ref. [57]) that is given by dui
dt ac
=X
j
mjαacv˜ijsig
ρij (ui−uj)rij ·∇iWij
|rij| , (3.55) where αac ∈[0,1] is the artificial conduction coefficient in analogy to a constant thermal conductivity. However, in contrast to Eq. (3.46) the internal energyu is used directly instead of the temperatureT.
3.1.8 Equations of state
There is one important key issue that has not been mentioned so far, although it has serious impact on the outcome of a simulation. That is the pressure calcula-tion and, if heat transfer is also considered, also the temperature. As discussed in Sec. 2.1.2, both are given by the so called equations of state (EOS), where
the most suitable EOS must be chosen individually as the case arises depending on the actual system under consideration. Incidentally, the sound speed is also derived from the EOS, since it can be calculated byp
∂P/∂ρ. However, even in review articles, a discussion of the topic is often missed out or, at least, only mini-mally discussed. Although, the importance for hydrodynamics is obvious and it is mainly controlling thermodynamic aspects. This is in particular true for complex fluids, such as in phase separating fluids. One possible reason may stem from the fact that for most astrophysical purposes the ideal gas EOS is a reasonably good approximation, for why alternatives may be considered as less important to discuss. Applying Eq. (2.5) to SPH yields
Pi= (γ−1)ui·ρi (3.56)
for a particle i, where γ = (f + 2)/f = CP/CV is the adiabatic index of the actual gaseous system. If the system under consideration is a liquid, pressure should increase at least with ρ2. Thus, the EOS should be replaced by another that better reproduces this dependency. One common EOS for the simulation of liquids is the so called Tait equation (see e.g. Ref. [58]),
Pi = ρ0ci
γ0
"
ρi
ρ0 γ0
−1
#
+P0, (3.57)
where ci is the sound speed. Pressure always pushes the system back to a be-forehand predicted normal pressure P0, which is due to the increase of pressure as the density deviates from the reference density ρ0. The strength is controlled by an exponent γ0, which is a numerical parameter and must not to be confused with the adiabatic index. Typical choices are γ0 ∈ {2,7}. However, this is only an empirically determined equation to fit experimental data of water, that is not based on a phenomenological approach, such as the ideal gas EOS. Moreover, it is just an isothermal EOS.
One EOS that is based on clear and simple physical assumptions (for a detailed discussion see Sec. 2.1.2), is the van der Waals equation of state (vdW-EOS). It can be understood to combine the linear ideal gas approach with a higher order density dependency as in the Tait equation. As it is also based on a phenomeno-logical approach is is an excellent choice for the simulation of liquid-vapor phase transitions. The vdW-EOS are described in detail in Sec. 2.1.2 and given by the two EOS (2.7) and (2.9). They read, applied to a particle i, as
Pi = ρikBTi
1−ρib−aρ2i, (3.58)
and
ui = kBTi
γ−1 −aρi, (3.59)
Chapter 3 Algorithms and methods
whereγis the adiabatic index andkB=kB/µmp, withkBthe Boltzmann constant, µ the mean molecular weight and mp the proton mass, and a, b are cohesive pressure and covolume, respectively, given in specific units (see Sec. 2.1.2).
As already mentioned, it is due to the second term in Eq. (3.58) that a neg-ative pressure can occur in a distinct density interval. From a classical point of view, there is no explanation for this behavior, e.g. this should never happen for an ideal gas. However, in the context of SPH this is a well known outstanding problem9, with a need to comply with artificial numerically controlled conditions.
Nevertheless, in the context of phase separation, especially regarding the dynamic evolution of phase transitions, it even seems to be intended in order to obtain an unstable region with negative pressures below the spinodal curve, just as in the vdW theory. Therefore, in this context there is no need to suppress negative pressures, but rather simply initiating the phase separation process.
9See Ref. [12] for a more detailed discussion on this, so called, tensile instability, which concerns the appearance of negative pressures in SPH. This must not be confused with the pairing instability, which concerns choosing the right smoothing length to prevent particle clumping, as discussed in Sec. 3.1.7.