Question about the algorithm - cfdemSolverPiso.C

Submitted by Tao on Fri, 02/17/2012 - 10:59

Hi all,
I am not quite understand about the algorithm implemented in the file cfdemSolverPiso.C. Here I am trying to explain my question clearly. In that file, we get the fluid velocity as : U = rUA*UEqn.H() at the beginning of the PISO loop. Generally speaking, the term H() represents the matrix coefficients of the neighbouring cells multiplied by their velocity and the unsteady term and all the sources except the pressure gradient. My question is that why we use this equation: U -= voidfraction*rUA*fvc::grad(p) - Ksl/rho*Us*rUA to calculate the exact velocity at the end of the PISO loop? I think we don't need the term Ksl/rho*Us*rUA since it has already been counted in H(). Am I right? Can anybody tell me why we implement the coupling algorithm as it is in cfdemSolverPiso.C?

Thanks.

Best Regards!

Tao

Bruno | Mon, 03/02/2015 - 21:48

The first step where U = rUA * UEqn.H() is a process similar to a Rhie-Chow interpolation to prevent chekkerboard effect with pressure or Us field. Therefore, this U field does not feel the presence of Us OR p. This is why div Us appears as a pressure source term.

Consequently, to recover the right velocity field, you need to carry out the final operation U -= voidfraction*rUA*fvc::grad(p) - Ksl/rho*Us*rUA. This ensures that second order accuracy in both space and time are recovered.
Cheers
BB