Riemann Particle Interaction Scheme

The Riemann problem can be defined as a category of initial value problems that involve a conservation equation and a piecewise data set with a single discontinuity.

In terms of SPH, you can solve a Riemann problem for each pair of particles, yielding a more isotropic particle distribution and smoother gradients in the solution. The following depiction is based on the work of C. Zhang. Given the Riemann problem definition, you start off by constructing a unit vector e i j between two particles i MathType@MTEF@5@5@+= feaahqart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVCI8FfYJH8YrFfeuY=Hhbbf9v8qqaqFr0xc9pk0xbb a9q8WqFfeaY=biLkVcLq=JHqpepeea0=as0Fb9pgeaYRXxe9vr0=vr 0=vqpWqaaeaabiGaciaacaqabeaadaqaaqaaaOqaaabaaaaaaaaape GaamyAaaaa@36F7@ and j MathType@MTEF@5@5@+= feaahqart1ev3aqatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVCI8FfYJH8YrFfeuY=Hhbbf9v8qqaqFr0xc9pk0xbb a9q8WqFfeaY=biLkVcLq=JHqpepeea0=as0Fb9pgeaYRXxe9vr0=vr 0=vqpWqaaeaabiGaciaacaqabeaadaqaaqaaaOqaaabaaaaaaaaape GaamOAaaaa@36F8@ . The two particles can then be assumed to be the left and right sides of the Riemann problem, such that the physical values of density ( ρ ), velocity ( U ) and pressure ( p ) are expressed as:

( ρ L , U L , p L ) = ( ρ i , v i ⋅ e i j , p i ) ( ρ R , U R , p R ) = ( ρ j , v j ⋅ e i j , p j )

Where, v is the velocity vector of the respective particle, which is illustrated in Figure 1. C. Zhang
Figure 1. Construction of the Riemann problem in SPH

The discontinuity is assumed to be located at the interface (Figure 1), half way between the particles. Assuming that particles are moving with arbitrary velocity vectors, the consequent solution of their motion is either a rarefaction or compression wave. However, there is a third wave solution at the interface denoted with * exponent, where U L ∗ = U R ∗ = U ∗ and p L ∗ = p R ∗ = p , as indicated by the dashed link in Figure 2.
Figure 2. Simplified Riemann Problem Solution in Time. This example shows three possible waves.C. Zhang


From these basic assumptions, you can derive the U * and p * values as:

U * = U a v g + 1 2 ( p L − p R ) p a v g C 0
p * = p a v g + 1 2 p a v g C 0 ( U L − U R )

Where, the a v g index denotes the mean value of the respective value between the two particles.

With these two solutions to the Riemann problem, you can modify the continuity equation and the pressure gradient term in the momentum equation as follows:

d p i d t = 2 p i ∑ j m j p j ( v i − v * ) ⋅ ∇ i W i j d v i d t = 2 ∑ j m j ( p * p i p j ) ∇ i W i j

Where, v * = U * e i j + ( v a v g − i j − U a v g e i j ) .

It is well-known that 1st order Riemann solvers such as the one illustrated here are dissipative. To that extent, nanoFluidX is using a limiter so that the Riemann solutions are provided only when the fluid is under compression. This effectively reduces the dissipation and has been shown to improve the results in various situations. C. Zhang
Note: Although the Riemann solver option is experimental, canonical test cases, simple geometries, and problems involving single phase runs such as tank sloshing or water management should still greatly benefit from the new formulation. However, industrial gearboxes and violent multiphase flows in complex moving geometries in general are likely to experience local instabilities. This will be addressed in the coming versions of nanoFluidX.
Density filtering is an integral part of the Riemann solver implementation. In SPH simulation, the consistency between mass, density, and occupied area cannot be strictly enforced. Several density filtering schemes have been developed to tackle this. In nanoFluidX, two versions are implemented: INCREMENTAL and INSTANT filtering.
INSTANT
The instant density filtering is basically a 0th-order Shepard filter. At each evaluation point (according to the number set by RM_freq_rho_init), the density is re-calculated by:
ρ i = ∑ j m j W i j ∑ j V j W i j MathType@MTEF@5@5@+= feaahqart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdi3aaS baaSqaaabaaaaaaaaapeGaamyAaaWdaeqaaOWdbiabg2da9maalaaa paqaa8qadaqfqaqabSWdaeaapeGaamOAaaqab0WdaeaapeGaeyyeIu oaaOGaamyBa8aadaWgaaWcbaWdbiaadQgaa8aabeaak8qacaWGxbWd amaaBaaaleaapeGaamyAaiaadQgaa8aabeaaaOqaa8qadaqfqaqabS WdaeaapeGaamOAaaqab0WdaeaapeGaeyyeIuoaaOGaamOva8aadaWg aaWcbaWdbiaadQgaa8aabeaak8qacaWGxbWdamaaBaaaleaapeGaam yAaiaadQgaa8aabeaaaaaaaa@4B60@
Where,
ρ i MathType@MTEF@5@5@+= feaahqart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdi3aaS baaSqaaiaadMgaaeqaaaaa@38CD@
Density of particle I
m j
Mass of particle j
V j
Volume of particle j
W i j
Weight calculated from the kernel function
INCREMENTAL
A modified 1st-order density filtering scheme is used as an additional option. The operator is calculated following:
ρ i = ∑ j m j W i j min 1 , ∑ j m j ρ * W i j MathType@MTEF@5@5@+= feaahqart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdi3aaS baaSqaaabaaaaaaaaapeGaamyAaaWdaeqaaOWdbiabg2da9maalaaa paqaa8qadaqfqaqabSWdaeaapeGaamOAaaqab0WdaeaapeGaeyyeIu oaaOGaamyBa8aadaWgaaWcbaWdbiaadQgaa8aabeaak8qacaWGxbWd amaaBaaaleaapeGaamyAaiaadQgaa8aabeaaaOqaa8qaciGGTbGaai yAaiaac6gadaqadaWdaeaapeGaaGymaiaacYcadaqfqaqabSWdaeaa peGaamOAaaqab0WdaeaapeGaeyyeIuoaaOWaaSaaa8aabaWdbiaad2 gapaWaaSbaaSqaa8qacaWGQbaapaqabaaakeaacqaHbpGCdaahaaWc beqaa8qacaGGQaaaaaaakiaadEfapaWaaSbaaSqaa8qacaWGPbGaam OAaaWdaeqaaaGcpeGaayjkaiaawMcaaaaaaaa@544A@
Where, ρ i MathType@MTEF@5@5@+= feaahqart1ev3aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr 4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9 vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdi3aaS baaSqaaiaadMgaaeqaaaaa@38CD@ is the density from solving the continuity equation before applying the filter.