Next: 8.2 Incompressible Flow with
Up: 8.1 Compressible Flow without
Previous: 8.1.2 Particle-in-Cell Method (PIC)
8.1.3 Smoothed Particle Hydrodynamics (SPH)
- PIC technique uses both Eulerian and Lagrangian elements.
Average density within a cell = number of point particles
in that cell.
- Can we represent the local fluid density without a grid?
Lucy [LUCY 77] and Gingold and Monaghan [GINGOLD 77,MONAGHAN 92]:
load each particle with a spatially extended
interpolation kernel
Average local density
= sum over the individual contributions.
Let
denote the interpolation kernel centered
around
. Then the local density at
is
|
(8.33) |
Generally, a property
is represented by
its ``smoothed particle estimate''
|
(8.34) |
Form of the interpolation kernel: Gaussian or polynomial
Example:
|
(8.35) |
with a width chosen such that the number of particles within
is in 2 dimensions and for 3 dimensions.
Now rewrite the Lagrangian equations of motion
8.8, 8.9 and 8.11.
in smoothed particle form.
Note:
In the momentum equation
,
interpolating and directly would not conserve
linear and angular momentum [MONAGHAN 92].
Use the identity
|
(8.36) |
and the SPH expressions for
and
to find
|
(8.37) |
with
.
If is Gaussian, this equation describes the motion of particle
under the influence of central pair forces
|
(8.38) |
The SPH equivalents of the other Lagrangian flow equations are
|
(8.39) |
where
, and
|
(8.40) |
Note: The density equation need not be integrated; just update all
positions , then invoke
equ. 8.33 to find
.
To update
the obvious relation
|
(8.41) |
might be used; a better way is
|
(8.42) |
with
. This relation
maintains angular and linear momentum conservation,
with the additional advantage that nearby particles will have similar
velocities [MONAGHAN 89].
To solve equs. 8.39, 8.37, 8.41 and
8.40, use any suitable algorithm (see Chapter 4).
Popular schemes are the leapfrog algorithm, predictor-corrector and
Runge-Kutta methods.
Example: Variant of the half-step technique ([MONAGHAN 89]):
Given all particle positions at time , the local density
at
is computed from
8.33. Writing equs. 8.37 and 8.40 as
|
(8.43) |
compute the predictors
|
(8.44) |
and
|
(8.45) |
Determine mid-point values of
,
and according to
|
(8.46) |
etc. From these, compute mid-point values of ,
and and insert these in correctors of the type
|
(8.47) |
Here is an overview of the method:
Smoothed particle hydrodynamics (SPH)
Note: SPH is not restricted to compressible inviscid flow.
Incompressibility may be introduced by using an equation of state
that keeps compressibility below a few percent [MONAGHAN 92],
and viscosity is added by an additional term in the
equations of motions for momentum and energy, equs. 8.37
and 8.40:
The artificial viscosity term is modeled in the following way:
|
(8.50) |
where is the speed of sound, is defined by
|
(8.51) |
and
and
etc.
This form of introduces the effects of shear and bulk viscosity.
The parameters and should be near
and for best results [MONAGHAN 92]. The quantity
prevents singularities for
. It should be
chosen such that
.
Additional features:
Thermal conduction may be included. See [MONAGHAN 89].
Interfaces: Introduce dummy particles on the far side of
the boundary. By picking the properties of these particles appropriately
one can mimick a free surface or a ``sticky'' solid boundary.
See Nugent and Posch [NUGENT 00] for
free surfaces, and Ivanov [IVANOV 00] for rough interfaces
Sample application:
``Rayleigh-Bénard'' convection
A fluid layer is heated carefully from below and cooled from above.
Formation of stable convective rolls transporting heat from
the bottom to the top. See the Figure for a match
between SPH and an Euler-type calculation [HOOVER 99]:
- Computing times comparable for both calculations
- Results are in good agreement
- Fluctuations in SPH (like in any particle-type calculation), none
in Euler
- SPH code is quite simple - similar to an MD program; Euler code
very massive
Comparison of Smoothed Particle Hydrodynamics with an
Eulerian finite-difference calculation. The density (above) and
temperature (below) contours for a stationary Rayleigh-Bénard
flow are shown. Left: SPH; right: Euler.
(From [HOOVER 99], with kind permission by the author)
Next: 8.2 Incompressible Flow with
Up: 8.1 Compressible Flow without
Previous: 8.1.2 Particle-in-Cell Method (PIC)
Franz J. Vesely Oct 2005
See also: "Computational Physics - An Introduction," Kluwer-Plenum 2001