[zx-81]

Hi,

I am working on fluid simulation code mainly intended for visualization/animation purposes and do not have a solid background in theoretical CFD so please excuse me if the following may be too basic/dumb kind of question:

For 'reasons of code legacy' I use collocated grids throughout my Naiver Stokes solver framework and I realize that this comes with stability problems because of the well known 'weak velocity-pressure coupling' problem (checkerboard effect)

I am aware of the common correction schemes (Rhie-Chow etc.) but I came up with the following Idea of my own that seams simpler to implement (at least in my case):

As the reason for the checker-board instability is the usage of simple central-differences for the velocity gradients used to calculate cell divergence ('right-hand-side') in the pressure correction step:

So why not use a slightly better gradient estimation technique instead: For example why not fit a quadratic spline over the 3x3x3 neighborhood at each cell and take the analytical spline derivative as the cell-centered gradient.
This way it is obvious that *each* grid point in the neighborhood, including the central one, contributes to the value of the gradient estimation.

So wouldn't this in effect eliminate the checkerboard-problem of the central-difference-based scheme where the middle point is completely ignored for the gradient ?

Furthermore this would be simple to implement on top of an existing solver because a 3x3x3 quadratic spline does not seem to introduce new 'boundary case problems' because we still only use a maximum local distance of one grid neighbor in each dircetion ...

It does also not seem to significatnly affect performance (although the overhead for gradient calculation increases substantiallly this is still a small fraction of the total run-time which is usually dominated by the pressure/poisson and advection steps.

On the other hand I did not found mentioned about this rather simple idea in connection with collocated grid correction schemes in the literature/web-search and there must be reason why this is not widely mentioned/used. So probably it is me with my limited CFD knowledge obviuously erring and/or missing something so that's why I ask in this forum what the experts think about this

Thanks in advance for any answers

zx-81

[FMDenaro]

Originally Posted by zx-81

Hi,

I am working on fluid simulation code mainly intended for visualization/animation purposes and do not have a solid background in theoretical CFD so please excuse me if the following may be too basic/dumb kind of question:

For 'reasons of code legacy' I use collocated grids throughout my Naiver Stokes solver framework and I realize that this comes with stability problems because of the well known 'weak velocity-pressure coupling' problem (checkerboard effect)

I am aware of the common correction schemes (Rhie-Chow etc.) but I came up with the following Idea of my own that seams simpler to implement (at least in my case):

As the reason for the checker-board instability is the usage of simple central-differences for the velocity gradients used to calculate cell divergence ('right-hand-side') in the pressure correction step:

So why not use a slightly better gradient estimation technique instead: For example why not fit a quadratic spline over the 3x3x3 neighborhood at each cell and take the analytical spline derivative as the cell-centered gradient.
This way it is obvious that *each* grid point in the neighborhood, including the central one, contributes to the value of the gradient estimation.

So wouldn't this in effect eliminate the checkerboard-problem of the central-difference-based scheme where the middle point is completely ignored for the gradient ?

Furthermore this would be simple to implement on top of an existing solver because a 3x3x3 quadratic spline does not seem to introduce new 'boundary case problems' because we still only use a maximum local distance of one grid neighbor in each dircetion ...

It does also not seem to significatnly affect performance (although the overhead for gradient calculation increases substantiallly this is still a small fraction of the total run-time which is usually dominated by the pressure/poisson and advection steps.

On the other hand I did not found mentioned about this rather simple idea in connection with collocated grid correction schemes in the literature/web-search and there must be reason why this is not widely mentioned/used. So probably it is me with my limited CFD knowledge obviuously erring and/or missing something so that's why I ask in this forum what the experts think about this

Thanks in advance for any answers

zx-81

if you want to discretize in such a way, then you must discretize also Div(Grad p) as same way as you do for Div V . That will produce a very large stencil for the pressure

[zx-81]

Originally Posted by FMDenaro

if you want to discretize in such a way, then you must discretize also Div(Grad p) as same way as you do for Div V . That will produce a very large stencil for the pressure

Thanks for the answer.

Yes I thought about this, but hoped to get away with just keeping the simple stencil for the pressure correction poisson equation. After all (at least in my projection based scheme) the pressure correction is only an approximation anyway and so I thought it might be sufficient to use my 27-point quadratic B-spline gradient discretization only for the gradients and only for the calculation of the right-hand-side of the pressure correction (and possibly for the pressure gradient but only for the result when adding it to the velocity after the poisson solve...
(somewhat similar to Chie-Row-Interpolation if I understood that correctly)

[zx-81]

After some more research I found a paper

"A velocity-pressure Navier-Stokes solver using a B-spline collocation method"
http://ctr.stanford.edu/ResBriefs99/botella.pdf

dealing with the same idea.

However, I'm afraid that the technical details of that paper are beyond my comprehension. In particular they mention that even with the B-spline discretization, the oscillation problem still occurs and they recommend a
kind of "staggered B-spline basis" construction to overcome it.

But this seems to be even more complicated than the regular staggered grid approch that I was trying to avoid wih the B-spline-approach ....

So the error must be in my original assumption that it is the central difference gradient discretization alone that is responsible for the oscillation instability because with B-spline gradients, we use not only two but 27 neighboring grid cells for each gradient value and so I thought that this should not be prone to checkerboard oscillations any more ... but what do I know

[FMDenaro]

but you want an exact or an approximate projection ?

[zx-81]

Originally Posted by FMDenaro

but you want an exact or an approximate projection ?

an approximate (but stable, non-oscillating) projection would be sufficient.

The simulation is only for visualization/animation purposes, but the grid resolution is high (up to 1000x1000x1000) and unfortunately, at resolutions higher than about 200^3 instabilities/checkerboard/striping artifacts occur unless I use agressive artificial damping of the velocity field (which of course negates the effect of the high resolution :-/

[FMDenaro]

are you sure this is not a symptom of numerical instability? maybe your time-step is to high for such a refined grid

[zx-81]

Originally Posted by FMDenaro

are you sure this is not a symptom of numerical instability? maybe your time-step is to high for such a refined grid

As I use an unconditionally stable semi-lagrangian advection scheme (cubic interpolation accuracy) it is normally not necessary to pay too much attention to the time step.
To the contrary: choosing too small a time step is even contraproductive with the semi-lagrangian schemes because it can lead to many repeated advection-interpolation steps 'within the same grid cell' and repeated interpolation = heavy blurring/numerical dissipation ...

So I know what I should really do is to switch the whole framework from
the 'collocated-grid/semi-lagrangian' to a 'staggered-grid/all-eulerian'
paradigm .... but this would mean to re-write almost all the code from scratch

However, on the other hand I really like the ability to take arbitrarily large time steps of the semi-lagrangian scheme (probably the only way to do a 1000^3 simulation on a normal workstation and the computational simplicity of the collocated grid ...

Concerning the co-located grid problem: maybe I schould take another look at the Rhie-Chow correction (but then again, wouldn't this mean to also change the discretization of the pressure gradient to a larger stencil in the poisson solver, which I try to avoid ?)

[FMDenaro]

I worked on this problem some time ago, if you want check details here are you. However, the cost is a second elliptic solver...

1. 
   A 3D second-order accurate projection-based Finite Volume code on non-staggered, non-uniform structured grids with continuity preserving properties: application to buoyancy-driven flowsInternational Journal for Numerical Methods in Fluids
   
   Volume 52, Issue 4, 10 October 2006, Pages: 393–432, F. M. Denaro
   Article first published online : 10 FEB 2006, DOI: 10.1002/fld.1185
   
2. 
   A non-diffusive, divergence-free, finite volume-based double projection method on non-staggered gridsInternational Journal for Numerical Methods in Fluids
   
   Volume 53, Issue 7, 10 March 2007, Pages: 1127–1172, A. Aprovitola and F. M. Denaro
   Article first published online : 27 SEP 2006, DOI: 10.1002/fld.1368
   

[zx-81]

Originally Posted by FMDenaro

I worked on this problem some time ago, if you want check details here are you.

Thank you very much ! this turned out to be really helpful.